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Abstract 

Currently,  there  are  four  Global  Navigation  Satellite  Systems  (GNSS)  either 
being  developed  or  in  existence- GPS,  GLONASS,  Compass,  and  Galileo.  Addition¬ 
ally,  there  are  several  Regional  Navigation  Satellite  Systems  (RNSS)  planned  or  in 
existence,  as  well  as  numerous  augmentation  systems  (which  require  a  GNSS  for  op¬ 
eration).  It  can  be  anticipated  that  there  will  be  interest  in  developing  additional 
independent  regional  navigation  satellite  systems  to  cover  areas  of  interest  to  partic¬ 
ular  countries  or  regions,  who  want  to  have  their  own  system. 

In  this  paper,  a  genetic  algorithm  is  used  in  an  effort  to  determine  near-optimal 
RNSS  constellations.  First,  a  cost  function  is  setup,  which  involves  a  weighted  com¬ 
bination  of  dilution  of  precision  (DOP)  values  and  percentage  availability  for  any 
number  of  receiver  locations  on  the  ground  (which  themselves  can  be  weighted).  Ef¬ 
fectively,  using  this  approach  it  is  easy  to  quantify  the  quality  of  coverage,  in  terms  of 
measurement  geometry,  over  a  specific  region  of  the  earth.  Next,  a  genetic  algorithm 
is  used  in  order  to  attempt  to  converge  to  the  lowest-cost  constellation  possible. 

Due  to  the  extremely  large  size  of  the  search  space,  there  appears  to  be  benefit 
to  constraining  the  search  space,  using  reasonable  knowledge  of  the  problem.  One 
example  is  to  constrain  the  constellation  to  include  two  geostationary  satellites.  This 
approach  was  taken,  and  slightly  better  costs  were  obtained  over  the  unconstrained 


IV 


case.  What  is  interesting  to  note,  however,  is  that  the  genetic  algorithm  tended 
to  result  in  generating  three  geostationary  satellites  (when  only  two  were  specified), 
strongly  suggesting  that  three  geostationary  satellites  is  the  correct  approach  for  the 
region  that  was  evaluated.  This  paper  shows  the  solutions  generated  by  genetic  algo¬ 
rithm,  and  conclusions  are  drawn  about  what  appear  to  be  the  best  constellations  for 
the  specified  region.  Overall,  the  genetic  algorithm  approach  appears  to  be  a  useful 
tool  for  constellation  optimization  for  RNSS  systems. 
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CONSTELLATION  DESIGN  OF  GEOSYNCHRONOUS 
NAVIGATION  SATELLITES  WHICH  MAXIMIZES 
AVAILABILITY  AND  ACCURACY  OVER  A  SPECIFIED 

REGION  OF  THE  EARTH 

I.  Introduction 

1.1  Satellite  Navigation  Background 

The  space  age  began  with  the  launch  of  first  artificial  satellite,  Sputnik  I  by  the 
Soviet  Union  in  1957.  The  launch  of  Sputnik  I  triggered  the  desire  to  discover  the 
details  of  the  satellite  orbit.  It  was  possible  to  determine  the  position  of  the  Sputnik 
by  the  radio  signals  transmitted  by  itself.  It  was  recognized  that  this  technique  could 
be  inverted:  if  the  satellite’s  orbit  was  already  known,  a  radio  receiver’s  unknown 
position  could  be  determined  accurately  by  using  the  same  type  of  measurements. 
Satellites  could  be  designed  specifically  for  the  purpose  of  providing  location  .  This 
idea  was  realized  in  the  Transit  Navigation  Satellite  System  by  the  U.S.  Navy.  Two 
virtually  identical  systems  were  developed  by  the  Soviet  Union:  Parus  for  the  Soviet 
Navy  and  Tsikada  for  merchant  ships  [23]. 

The  success  of  Transit  motivated  the  U.S.  Navy  and  the  Air  Force  for  follow-on 
programs  for  space-based  navigation  systems  which  resulted  in  the  Global  Positioning 
System  (GPS).  The  Soviet  Union  undertook  to  develop  a  similar  system  called  GLObal 
NAvigation  Satellite  System  (GLONASS).  The  GPS  is  fully  operational  and  has  var- 
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ious  application  areas  today.  After  the  dissolution  of  the  Soviet  Union,  the  Russian 
Federation  who  has  assumed  the  responsibility  for  GLONASS  was  not  able  keep  the 
GLONASS  fully  functional  due  to  economical  problems.  The  GLONASS  needs  more 
satellites  to  compete  with  the  GPS.  Therefore,  being  provided  the  required  budget, 
GLONASS  is  planned  to  have  global  coverage  at  the  end  of  2009  [24,28]. 

Another  project  in  development  is  the  European  system,  Galileo,  which  has 
launched  its  first  satellite  in  2005.  The  full  constellation  was  planned  to  be  operational 
by  2009.  However,  it  is  unlikely  that  full  operational  capability  will  be  available  until 
2012-2013  [27], 

There  are  also  Satellite  Based  Augmentation  Systems  (SBASs)  which  are  es¬ 
sentially  extra  satellites  transmitting  signals  intended  to  address  shortcomings  of 
GNSS  for  enhanced  accuracy,  availability,  and  integrity  for  civil  aviation  users.  The 
GPS  augmentation  system,  Wide  Area  Augmentation  System  (WAAS),  performs  well 
in  the  US  and  it  is  widely  used.  The  European  augmentation  system,  European 
Geostationary  Navigation  Overlay  System  (EGNOS)  (also  includes  corrections  for 
GLONASS),  is  becoming  mature  and  will  be  ready  for  certification  for  Safety  of 
Life  operations  in  mid  2008  over  Europe.  Japan  prepares  for  operations  with  the 
augmentation  system  MTSAT  (Multi-Functional  Transport  SATellite)  based  Satellite 
Augmentation  System  (MSAS),  interoperable  with  EGNOS  and  WAAS.  Japan  is  pro¬ 
ceeding  with  the  regional  GPS  augmentation,  Quasi-Zenith  Satellite  System  (QZSS), 
which  will  have  the  first  satellite  launched  in  2009.  India  is  proceeding  with  its  re- 
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gional  augmentation  system  which  is  called  GPS  and  GEO  Augmented  Navigation 
(GAGAN)  [28], 

Regional  Navigation  Satellite  Systems  (RNSSs)  are  intended  to  provide  naviga¬ 
tion  service  over  a  nation  or  region.  Their  satellite  constellations  are  much  smaller 
than  GNSSs,  perhaps  only  5-7  satellites  (depending  on  the  specific  orbit  configura¬ 
tion).  They  may  be  considered  as  first  steps  to  GNSSs  [27]. 

Japan  has  intentions  to  build  a  regional  system,  Japanese  Regional  Naviga¬ 
tion  Satellite  System  (JRANS),  which  will  consist  of  3  QZSS  satellites,  the  MTSAT 
satellite  and  three  additional  satellites  in  high  earth  orbit  [27]. 

China  has  developed  an  RNSS  which  is  called  BeiDou-1  and  it  is  reported  to 
have  an  accuracy  of  dozens  of  meters  [24] .  China  participates  in  the  Galileo  program 
and  is  developing  a  global  satellite  navigation  system  BeiDou-2  (or  COMPASS )  as  a 
follow  up  of  the  regional  BeiDou-1  [28]. 

The  Indian  government  has  started  a  project  of  an  independent  regional  naviga¬ 
tion  satellite  system  which  is  planned  to  be  completed  by  2011.  The  Indian  Regional 
Navigation  Satellite  System  (IRNSS)  will  have  a  seven  satellite  constellation  [15]. 

There  are  also  plans  for  completely  commercial  satellite  navigation  systems. 
Some  companies,  like  Lockheed  Martin  in  cooperation  with  NEC  of  Japan  and  Alcatel 
of  France  and  Boeing,  have  plans  for  their  own  augmentation  systems.  It  is  appealing 
to  own  a  satellite  based  navigation  system  due  to  the  political  and  economic  benefits. 
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Just  as  satellite  based  communications  was  commercialized,  governments,  companies 
and  multinational  conglomerates  may  be  competing  to  exploit  those  benefits  [18]. 

Global  Navigation  Satellite  System  (GNSS)  is  the  generic  name  for  the  class 
of  systems  of  which  GPS  is  currently  the  only  working  system.  The  technology  has 
advanced  significantly  and  the  world  political  order  has  shifted  since  GPS  was  de¬ 
signed  in  the  1970s.  The  U.S.  wants  to  maintain  leadership  among  GNNSs  and  has 
plans  to  develop  GPS  Ill  which  is  a  more  capable  system.  The  GLONASS  is  planned 
to  be  revived,  and  the  Galileo  has  the  momentum  and  the  funds  to  be  operational. 
The  GNSS  world  is  going  towards  more  competition  and  cooperation.  It  is  advanta¬ 
geous  to  have  more  satellites  and  better  signal  quality  of  multiple  systems  for  users. 
Major  components  are  the  USs  modernized  GPS  and  planned  GPS-III,  the  revital¬ 
ized  GLONASS,  and  Europes  planned  GALILEO  system.  In  addition,  the  SBASs 
and  RNSSs  will  provide  extra  satellites  and  signals  to  the  GNSS  ’’mix”.  The  user 
community  will  benefit  more  from  modernized  GPS,  GLONASS,  Galileo,  SBASs  and 
RNSSs  [24,27], 

The  future  of  GNSS/RNSS  is  promising.  More  satellites  and  more  signals  will 
be  welcome  by  many  users.  Better  performance  can  be  achieved  by  integrating  all 
these  satellite  systems  into  one  Global  Navigation  Satellite  System  of  Systems.  The 
benefits  of  satellite-based  positioning,  navigation  and  timing  offered  by  compatibility 
and  interoperability  of  all  GNSSs,  SBASs  and  RNSSs  are  far  greater  than  what  any 
individual  system  can  provide. 


4 


1.2  Problem  Statement 


As  stated  in  Section  1.1,  having  a  satellite  navigation  system  is  politically  and 
economically  beneficial.  On  the  other  hand,  building  and  maintaining  a  satellite 
navigation  system  is  extraordinarily  expensive  depending  on  the  size  of  the  satellite 
constellation.  In  order  to  minimize  the  cost  and  maximize  the  benefit  of  such  a 
system,  the  requirements  for  the  coverage  area  and  the  expected  performance  of  the 
system  should  be  determined.  Based  on  those  requirements,  an  optimized  satellite 
constellation  that  meets  the  requirements  can  be  designed.  There  has  been  some 
studies  to  provide  solutions  for  both  regional  and  global  coverage  requirements. 

The  previous  research  on  regional  and  global  navigation  satellite  systems  has 
proposed  solutions  which  are  exploiting  the  combinations  of  common  orbits,  like  Low 
Earth  Orbits  (LEO),  Geostationary  Orbits  (GEO),  Geosynchronous  Orbits  (GSO), 
and  Medium  Earth  Orbits  (MEO),  as  solutions  [7,8,17,20,26]. 

The  previous  research  on  the  satellite  constellation  design  using  Genetic  Al¬ 
gorithms  (GA)  has  showed  the  power  of  the  GA  in  satellite  constellation  optimiza¬ 
tion  [8,10-12,30-32],  Previous  research  used  GA  to  optimize  combinations  of  the 
above  mentioned  orbits  considering  the  distribution  of  the  orbital  planes  on  the  equa¬ 
tor  and  the  distribution  of  the  satellites  on  those  orbits.  While  much  work  has  been 
done  on  constellation  optimization,  very  little  has  been  done  to  design  a  constellation 
for  navigation  purposes  to  serve  within  a  specified  region  of  the  world.  The  goal  of 
this  research  is  to  design  a  constellation  of  geosynchronous  navigation  satellites  which 
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maximizes  availability  and  accuracy  over  a  specified  region  of  the  earth  by  using  the 
GA  for  optimization.  In  order  to  accomplish  this  goal;  a  simulation  tool  is  developed 
to  evaluate  the  navigation  performance  of  satellite  constellations  in  which  the  satel¬ 
lites  have  GEO  period  as  a  constant  and  the  other  orbital  parameters  as  variables, 
a  genetic  algorithm  is  designed  to  find  an  optimal  or  a  near-optimal  solution  to  the 
problem.  The  objective  is  to  find  a  solution  which  provides  100%  availability  with 
the  best  available  accuracy. 

1.3  Scope  and  Assumptions 

Simplifying  assumptions  are  needed  to  focus  the  research  in  a  direction  so  that 
conclusions  can  be  reached  without  making  the  simulations  overly  complex.  Several 
assumptions  were  made  to  limit  the  scope  of  the  research: 

•  The  simulation  of  the  satellite  orbits  for  optimization  is  based  on  two-body 
equations  of  motion.  For  the  analysis  of  the  final  candidate  solutions,  J2  per¬ 
turbation  is  included  to  make  it  more  realistic. 

•  The  performance  of  the  simulated  navigation  system  constellations  will  be  eval¬ 
uated  considering  the  accuracy  of  the  position  estimate  and  availability  of  a 
given  level  of  accuracy  . 

•  The  measurement  accuracy  assumed  to  be  consistent  for  a  given  system  (a  rea¬ 
sonable  assumption  commonly  made  with  satellite  navigation  systems).  There¬ 
fore,  the  accuracy  of  the  position  estimate  is  assumed  to  depend  only  on  the 
geometric  quality  of  the  constellation. 
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•  The  candidate  region  for  the  navigation  service  is  as  depicted  in  Figure  1.1. 


Figure  1.1:  The  candidate  region  for  the  navigation  service 


1.4  Related  Research 

Romay-Merino,  Cobo,  and  Herraiz-Monseco  studied  the  design  of  high  perfor¬ 
mance  and  cost  efficient  constellations  for  future  global  navigation  systems.  They 
analyzed  three  types  of  constellations:  GSO/GEO  constellations,  GPS  like  constella¬ 
tions,  and  hybrid  constellations  (GPS/GSO/GEO).  They  used  an  adaptive  random 
search  method  to  optimize  these  constellations  in  order  to  minimize  the  worst  perfor¬ 
mance  values  even  in  case  of  satellite  failures.  They  used  a  geographical  distribution 
of  weighting  function  to  have  the  constellation  provide  better  performance  over  highly 
populated  areas.  The  constraints  on  satellite  constellations  that  they  applied  include 
all  satellites  at  geosynchronous  altitude,  all  satellites  in  inclined  planes  have  the  same 
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inclination,  fixed  number  of  inclined  planes,  several  satellites  in  the  same  ground  track 
and  fixed  number  of  GEO  and  GSO  orbits. 

The  results  of  optimization  of  GSO/GEO  constellation  indicated  that  at  least 
35  satellites  are  required  to  accomplish  the  requirements.  They  also  reported  that 
failure  of  a  satellite  in  this  kind  of  a  constellation  would  significantly  deteriorate 
the  performance  of  navigation  service  over  the  area  covered  by  the  failed  satellite. 
They  also  performed  an  optimization  of  a  regional  GSO/GEO  constellation  to  provide 
better  performance  over  Europe.  They  found  that  at  least  22  satellites  are  required 
to  satisfy  the  requirements  over  Europe.  They  made  another  optimization  for  GPS- 
like  constellations  and  concluded  that  at  least  45  satellites  are  required  to  satisfy 
the  requirements  and  the  performance  is  about  the  level  of  performance  provided 
by  the  35  satellite  GSO/GEO  constellation  [8].  Finally,  they  optimized  the  hybrid 
(GPS/GSO/GEO)  constellations  and  concluded  that  at  least  42  satellites  are  required. 
The  overall  conclusions  for  their  study  as  follows  [8]: 

•  The  use  of  an  optimization  tool  allows  to  design  cost  efficient  and  high  perfor¬ 
mance  navigation  satellite  constellations. 

•  GSO/GEO  constellations  can  be  designed  to  provide  better  performance  over 
desired  regions  by  using  geographically  distributed  weighting  functions. 

•  A  regional  satellite  navigation  constellation  using  GSO/GEO  can  easily  be 
transformed  into  a  global  constellation  easily. 


Renault  discussed  the  properties  of  regional  navigation  constellations  which  are 
combinations  of  GEO  and  GSO  orbits.  He  stressed  that  synchronous  orbits  are  the 
best  suited  orbits  for  a  regional  service,  since  their  visibility  over  a  regional  area  is 
maximum.  Therefore  he  suggested  that  a  regional  constellation  should  be  composed 
of  both  GSO  and  GEO  satellites.  He  analyzed  the  performance  of  both  GEO  and 
GSO  based  constellations.  He  noted  that  GEO  based  constellations  had  better  ver¬ 
tical  performance  and  GSO  based  constellations  had  better  horizontal  performance. 
He  concluded  that  three  GEOs  and  three  GSOs  on  three  planes  define  a  minimal 
navigation  constellation  for  a  mid-latitude  region  [26]. 

Micheau  and  Thiebolt  studied  the  satellite  constellation  design  for  navigation 
needs.  They  specifically  studied  LEO  constellations  with  complementary  GEO  satel¬ 
lites.  The  low  cost  properties  of  LEO  satellites  made  them  appear  to  be  a  sensible 
choice.  They  emphasized  maximizing  the  position  accuracy  over  Europe.  They  used 
both  analytical  methods  and  genetic  algorithms  to  find  an  optimal  constellation  so¬ 
lution  to  meet  the  civil  aviation  specifications  over  Europe.  The  comparison  of  the 
results  obtained  by  different  methods  showed  that  the  the  results  of  Genetic  algo¬ 
rithms  were  the  only  one  that  met  all  requirements  [20]. 

Carnebianca  studied  regional  to  global  satellite  based  navigation  systems.  He 
suggested  to  create  a  constellation  with  two  GEO  satellites  and  two  or  three  satellites 
on  inclined  orbit  of  12  or  24  hour  period  for  a  regional  navigation  service.  The  three 
major  trade-off  elements  considered  in  his  study  are  number  of  GEO  satellites,  their 
location  and  equatorial  spacing  and  the  properties  of  inclined  satellites  such  as  number 
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of  inclined  satellites,  their  apogee  loop  location  and  their  phasing  relationship.  He 
analyzed  different  combinations  of  GEO,  Highly  Elliptical  Orbit  (HEO)1  and  Tundra2 
orbits  to  create  a  global  navigation  constellation.  He  emphasized  the  advantage  of 
GEO  satellites  for  serving  a  region  of  the  earth.  He  stated  that  the  use  of  GEO 
satellites  allows  a  regional  service  with  the  potential  to  be  expanded  gradually  with 
additional  satellites  [7]. 

Frayssinhes  investigated  satellite  constellation  geometries  with  genetic  algo¬ 
rithms  for  data  collection  purposes.  He  stressed  the  efficiency  of  genetic  algorithms  in 
optimizing  satellite  constellations  and  noted  the  effectiveness  of  the  GA  for  optimiz¬ 
ing  add-on  satellites  to  the  GPS.  He  proposed  relying  on  real  encoding  instead  of  a 
binary  encoding  mechanism.  He  used  the  GA  to  distribute  the  orbital  planes  on  the 
equator  and  to  distribute  the  satellites  on  the  orbital  planes  (i.e. ,  optimization  of  the 
right  ascension  of  ascending  node  and  the  initial  mean  anomaly).  He  concluded  that 
the  GA  with  a  binary  encoding  mechanism  was  not  suitable  to  manage  circular  con¬ 
figurations  that  are  both  distribution  of  satellites  in  an  orbital  plane  and  distribution 
of  the  orbit  ascending  nodes  over  the  equator  [12], 

Ferringer  and  Spencer  studied  satellite  constellation  design  optimization  via 
multiple-objective  evolutionary  algorithms  (Non-Dominated  Sorting  Genetic  Algo¬ 
rithm  II)  for  remote  sensing  purposes.  They  made  optimization  of  a  six-satellite 

1HEO:  period  =  12  hrs,  apogee/perigee  altitude  =  39105/1250  km,  inclination  =  63.45 

2Tundra:  period  =  24  hrs,  apogee/perigee  altitude  =  46340/25231  km,  inclination  =  63.45 
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constellation  considering  the  distribution  of  orbits  on  the  equator,  distribution  of 
satellites  on  those  orbits,  and  the  same  inclination  value  for  each  orbit  [11], 

Ely,  Crossley  and  Williams  introduced  a  method  for  applying  a  GA  approach  for 
multi-objective  design  of  constellations  for  zonal  coverage.  Their  study  has  shown  that 
genetic  algorithms  and  eccentric  orbits  are  useful  in  designing  efficient  constellations. 
Their  approach  provides  an  opportunity  to  identify  constellation  designs  with  elliptic 
orbits  that  require  less  launch  effort  and  cost  less  than  designs  using  circular  orbits 

[30]. 

1.5  Thesis  Overview 

In  Chapter  two,  a  background  of  position  estimation  in  GPS  and  Genetic  Al¬ 
gorithms  is  presented.  In  Chapter  three,  the  structure  of  the  simulation,  the  imple¬ 
mentation  of  the  GA,  and  the  orbit  propagation  are  described.  In  Chapter  four,  the 
data  generated  by  the  simulation  is  analyzed.  In  the  chapter  Eve,  the  results  of  the 
simulation  and  the  possible  benefits  the  study  are  discussed. 
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II.  Background 


2. 1  Overview 

The  purpose  of  this  chapter  is  to  acquaint  the  reader  with  the  key  concepts  used 
in  this  research.  In  Section  2.2,  a  brief  explanation  of  GPS  and  its  segments  is  covered. 
More  information  about  the  GPS  can  be  found  in  [24].  In  Section  2.3,  definition  of 
common  terms  to  describe  the  performance  of  a  satellite  navigation  system  is  given. 
In  Section  2.4,  the  position  estimation  in  GPS  and  the  Dilution  of  Precision  concept 
are  explored.  In  Section  2.5,  a  basic  definition  of  GA  and  common  operators  are 
provided. 

2.2  Global  Positioning  System 

GPS,  declared  operational  in  1995,  consists  of  three  segments:  the  user  segment, 
the  control  segment,  and  the  space  segment  [24],  The  user  segment  consists  of  all  kinds 
of  civilian  and  military  receivers.  The  control  segment  deals  with  the  management  of 
satellite  operations.  The  space  segment  comprises  the  satellites.  GPS  space  segment 
nominally  consists  of  24  satellites  in  nearly  circular  orbit  with  a  radius  of  26560  km, 
a  period  of  approximately  12  hours,  and  stationary  ground  tracks.  The  satellites 
are  arranged  in  six  equally  spaced  orbital  planes  with  an  inclination  of  55  degrees. 
Each  plane  hosts  four  satellites  unevenly  distributed  around  the  plane.  This  nominal 
constellation  assures  that  a  receiver  anywhere  in  the  world  would  see  at  least  four 
satellites  at  all  times.  Figure  2.1  illustrates  the  nominal  GPS  constellation.  In  order 
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to  provide  better  performance  to  the  user  community,  currently  31  functional  GPS 


satellites  are  in  space  (as  of  January  16,  2008)  [2], 


Figure  2.1:  Nominal  GPS  Constellation  [3]. 

2.3  Performance  of  a  Navigation  System 

The  performance  evaluation  of  a  navigation  system  is  usually  described  in  terms 
of  accuracy,  availability,  and  integrity. 

Accuracy  is  defined  as  the  difference  between  the  calculated  position  at  any 
given  time  and  the  true  position.  The  accuracy  of  a  navigation  system  depends  on 
the  accuracy  of  the  navigation  signal  received  by  the  user  and  the  relative  geometry 
of  the  navigation  satellites  of  which  signals  are  being  received  [26]. 

Availability  is  the  ability  of  a  navigation  system  to  provide  position  information 
whenever  it  is  needed  by  the  user.  The  accuracy  of  a  navigation  system  over  the  time 
is  not  constant.  When  the  accuracy  of  the  system  is  better  than  a  specified  value,  the 
system  is  called  available  [26]. 
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Integrity  is  the  ability  of  a  navigation  system  to  provide  timely  warnings  to 
users  or  to  shut  itself  down  when  the  accuracy  is  worse  than  a  specified  value  [26]. 

2.f  Position  Estimation  and  Dilution  of  Precision  (DOP) 

A  GPS  receiver  needs  signals  from  at  least  four  satellites  to  estimate  its  position. 
The  GPS  satellites  transmit  data  signals  which  includes  the  information  of  their  own 
position  in  space  and  the  time  of  transmission.  The  GPS  receiver  calculates  its  range 
to  each  visible  satellite  by  using  time-of- arrival  (TOA)  method.  In  order  to  utilize  the 
TOA,  the  GPS  receiver  must  have  the  information  of  the  time  of  signal  transmission, 
the  speed  of  signal  propagation,  and  the  time  of  signal  reception.  The  speed  of  signal 
propagation  for  GPS  signals  is  the  speed  of  the  light,  which  is  already  known.  The 
time  of  the  signal  transmission  is  provided  in  the  GPS  signal.  Although  satellites  have 
very  accurate  atomic  clocks,  the  time  information  is  not  as  precise  as  GPS  receiver 
needs.  However,  the  correction  for  satellite  clock  error  is  also  provided  in  the  GPS 
data  signal.  The  time  of  signal  reception  is  determined  by  the  receiver  clocks  which 
are  generally  inexpensive  quartz  clocks  which  has  relatively  less  stability  than  GPS 
satellite  clocks.  The  receiver  clock  error  is  not  a  problem,  because  it  can  be  explicitly 
solved  for  by  the  fourth  satellite  signal.  Therefore,  a  GPS  receiver  needs  at  least  four 
satellite  signals  to  determine  its  own  position. 

Figure  2.2  illustrates  the  vectors  from  the  receiver  to  the  satellites  in  Cartesian 
coordinate  system,  where  ( xu,yu,zu )  is  the  receiver  position  vector  and  (. Xk,yk,Zk )  are 
the  position  vectors  of  kth  satellite  for  k  —  1,  2,  ..,  n. 
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Figure  2.2:  Receiver-to-satellite  vectors 
In  order  to  determine  receiver  position  in  three  dimensions  ( xu,yu,zu )  and  the 
clock  offset  tu,  pseudorange  measurements  should  be  made  to  four  satellites  which 
can  be  written  as  follows  (assuming  the  only  error  is  receiver  clock  error): 

Pk  =  ||xfe  -  xJI  +  ctu  (2.1) 

where  pk  is  the  pseudorange  between  the  receiver  and  the  kth  satellite,  x*.  is  the  kth 
satellite  position  vector,  and  xu  is  the  receiver  position  vector.  Equation  2.1  can  be 
expanded  as: 


Pk 


=  V (xk  ~  xu)  +  (yk  ~  yu)  +  {zk  ~  Zu)  +  ctv 


f  Vut  Zui  tu) 


(2.2) 
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These  nonlinear  equations  can  be  solved  by  employing  iterative  techniques  based 
on  linearization.  If  the  approximate  receiver  position  is  known,  the  offset  of  the  true 
position  (xu,  yu ,  zu)  can  be  denoted  from  the  approximate  receiver  position  (xu,  yu,  zu) 
by  a  displacement  (Axu,  Ayu,  Azu) .  The  position  offset  (Axu,  Ayu,  Azu)  can  be  ob¬ 
tained  as  linear  functions  of  the  known  coordinates  and  pseudorange  measurements 
by  expanding  Equation  2.2.  By  using  the  approximate  position  location  (xu,yu,zu) 
and  time  bias  estimate  tu,  an  approximate  pseudorange  can  be  calculated: 


Pk 


=  \/ (xk  ~  xu)  +  (yk  -  yu)  +  (zk  *  zu)  +  civ 


f  Vui  ^Ul^u) 


(2.3) 


The  receiver  position  and  receiver  clock  offset  is  considered  to  consist  of  an 
approximate  and  an  incremental  component: 


xu 


Vu 

Zu 


xu  +  A  xu 
Vu  A  yu 
A  zu 


t 


u 


tu  T  Atu 


Therefore,  it  can  be  written  as: 


(2.4) 
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f  Vl Li  ZlL")  tu)  f  (%U  “1“  Vu  “I-  ^Vw)  %U  ^ZlLl  iu  ^tu) 


(2.5) 


Equation  2.5  can  be  expanded  about  the  approximate  point  and  associate  pre¬ 
dicted  receiver  clock  offset  (xu,yu,  zu,ctu)  using  a  Taylor  series  expansion: 


f  {xu  T  ^xu,  Vu  T  zu  T  ^Zui  iu  T  Atu)  f  (^u)  Vui  z ui  tu) T 


(2.6) 


9/(  Urn  in, K) 9/( £»,  ft.,  z«,  < J  .  9/( i», y»,  5», t») ^ _ 


OZq, 


df(xu,yu,zu,tu ) 

- x - Aru  +  higher  order  terms 

dtu 

In  order  to  eliminate  nonlinear  terms,  the  higher  order  terms  are  truncated.  The 
partial  derivatives  evaluate  as  follows: 


df{xu,yu,zu,tu ) 
dxu 

df(xu,yu,zu,iu) 

dyu 

df(xu,yu,zu,iu) 

dzu 

df(xu,yu,zu,iu) 

dtu 


Xk  Xu 

fk 

Vk  Vu 

fk 

Zk  ~  Zu 
fk 


=  c 


(2.7) 


17 


where  fk  is  the  magnitude  of  the  vector  from  the  approximate  user  position  to  the  kth 
satellite  and  defined  as: 


fk  =  \J ( xk  -  xu)2  +  (yk  -  yu )2  +  (zk  -  zu )2  (2.8) 

Substituting  Equations  2.4  and  2.7  into  Equation  2.6  yields: 


Pk 


pk 


-Ax, 


rk 


Vk  ~  Vu 

rk 


Ayu  ~ 


zk  ~  zu  „ 

— - - A  zu  +  cAtu 

rk 


(2.9) 


Equation  2.9  is  the  linearized  form  of  Equation  2.3  with  respect  to  the  unknowns 
Axu,  Ayu,  A zu,  and  A tu.  Rearranging  the  the  expression  with  the  known  quantities 
on  the  left  yields: 


A  .  &K  tfU  A  . 

Pk-  pk  =  7 - Axu  H - 7 - A  yu  + 

rk  rk 


rk 


-A zu  -  cAtu 


12.10) 


By  introducing  new  variables,  Equation  2.10  can  be  simplified  as: 


A  p  = 

Pk  Pk 

& xk 

Xk  xu 

fk 

dyk 

yk  -  i/u 

fk 

%k 

fk 

(2.11) 
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In  these  equations,  (axk,  dyk,  azk)  defines  the  unit  vector  pointing  from  the  approxi¬ 
mate  user  position  to  the  kth  satellite.  Substituting  the  new  variables  into  Equation 
2.10  yields: 


A Pk  ^xkAxu  T  &ykAyu  T  &zkAzu  cAtu 

The  four  unknowns  Axu ,  A yu,  A zu,  and  A tu  can  be  solved  for  by  making  ranging  mea¬ 
surements  to  at  least  four  satellites,  and  then  solving  the  following  linear  equations: 


Api  =  axlAxu  +  aylAyu  +  azlAzu  -  cAtu  (2.12) 

Ap2  aX'2Axu  T  ciy2Ayu  T  aZ2Azu  cAtu 
A p3  =  ax3Axu  +  dy3Ayu  +  dz3Azu  -  cAtu 


A pn  dxnAxu  T  aynA yu  T  dznAzu  cAt 


These  equations  can  be  expressed  in  matrix  form  as: 


- 

" 

< 

C^xl  dyl  zl  1 

Ap2 

ax  2  ay2  az2  —  1 

Ap3 

H  = 

& xS  & y3  ^z3  1 

> 

3 

(Lxn  CLyn  zn  1 

Ax  = 


Axu 
A  yu 
A  zu 

cA  tu 
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where  A p  is  pseudorange  difference  vector,  H  is  measurement  matrix  and  Ax  is  user 


position  displacement  vector.  This  matrix  equation  can  be  written  as: 


A  p  =  H  Ax 


(2.13) 


where  there  are  three  cases: 

1.  n  <  4  :  Underdeter  mined  case 

•  Cannot  solve  for  Ax 

2.  n  —  4  :  Uniquely  determined  case 

•  One  valid  solution  for  Ax 

•  Solved  by  calculating  HU1  (Ax  =  H_1Ap) 

3.  n  >  4  :  Overdeter  mined  case 

•  No  perfect  solution  that  solves  equation 

•  Least-squares  techniques  can  be  used  to  make  an  estimate 

The  solution  for  the  third  case,  which  is  more  common,  based  on  the  least-squares 
solution  is  as  follows: 


Ax  =  (HrH)-1HTAp 


(2.14) 


Once  the  Ax  is  calculated,  the  receiver’s  coordinates  xu,  yu,  zu  and  the  receiver 
clock  offset  tu  are  calculated  using  Equation  2.4.  This  linearization  works  well  as 
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long  as  the  displacement  (Axu,  Ayu,  Azu)  is  within  the  close  proximity  of  the  known 
approximate  receiver  position.  When  the  position  estimation  error  is  above  an  ac¬ 
ceptable  value,  taking  the  solution  as  new  approximate  known  receiver  position  and 
repeating  the  linearization  will  provide  a  better  estimate.  The  iteration  process  may 
go  on  until  the  desired  accuracy  is  reached  [16]. 

The  position  estimation  accuracy  also  depends  on  the  relative  geometry  of  the 
location  of  the  receiver  and  the  satellites  which  is  represented  by  H  matrix.  In  the 
beginning  of  the  position  estimation  calculation,  the  only  error  was  assumed  to  be  the 
receiver  clock  error.  In  a  real  world  situation,  the  receiver-to-satellite  measurements 
are  corrupted  by  independent  errors,  such  as  measurement  noise,  satellite  position 
error  (ephemeris  error),  multipath  error.  In  other  words,  A p  includes  all  of  the 
measurement  errors,  and  these  errors  are  amplified  by  (H7  H)-1  matrix  when  solving 
for  Ax.  Therefore,  a  badly  distributed  satellite  geometry  reduces  the  accuracy  of  the 
position  estimation  which  is  called  Dilution  of  Precision  (DOP).  The  good  and  poor 
DOP  situation  is  illustrated  in  Figure  2.3. 

The  DOP  concept  is  a  way  to  characterize  the  impact  of  measurement  geometry 
on  the  position  solution.  By  defining  the  covariance  matrix  of  measurements  (Cp)  and 
the  covariance  matrix  of  calculated  position  and  clock  error  (Cx)  as: 
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(a) 


(b) 

Figure  2.3:  (a)  Well  distribution  of  satellites  provides  a  good  DOP. 

(b)  Clustered  satellites  give  poor  DOP  values. 


C„  = 


c,  = 


< 

°rpiP2 

^PiPn 

<JP2Pl 

ap2 

<JP2Pn 

°PnPl 

&PnP2 

ax 

& xuyu 

@xuzu 

®xu  5tu 

^VuXu 

< 

<7yu  Z'u, 

^VuStu 

®zuxu 

OzuTJu 

cr \ 

&zu8tu 

^  5tuXu 

aStuyu 

&5tuzu 

a5tu 

(2.15) 


(2.16) 


and  making  the  following  assumptions; 


•  Measurement  errors  are  zero-mean 

•  Measurement  errors  have  a  Gaussian  distribution 
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•  All  measurements  have  the  same  variance 


•  Measurement  errors  are  uncorrelated 


the  following  definitions  can  be  made: 


c„  =  V 


(2.17) 


Cx  -  i  H  'Hi  V 


(2.18) 


The  matrix  (H7  H)  1  is  called  the  DOP  matrix  and  directly  relates  the  measurement 


errors  to  position  errors.  In  Equation  2.18  the  H  is  defined  in  Earth-Centered  Earth- 
Fixed  (ECEF)  frame,  but  DOP  values  should  describe  errors  in  geodetic  (local-level) 
coordinate  frame.  Therefore,  the  H  matrix  should  be  modified  so  that  the  unit  LOS 
vectors  are  expressed  in  the  local-level  frame.  When  the  H  matrix  is  thus  modified 
as  H(;,  the  Equation  2.18  is  redefined  as: 


cx  =  y;- 


(2.19) 


where 
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O en 

O eu 

®ebtu 

&ne 

Ol 

®nu 

&n5tu 

&ue 

®un 

ol 

&u5tu 

O Stue 

&5tun 

&StulL 

astu 

(2.20) 


It  is  desirable  to  characterize  the  Cx  matrix  using  a  number.  The  Geometric 
Dilution  of  Precision  (GDOP)  can  be  calculated  directly  from  the  modified  DOP 
matrix; 


Du 

D\2 

D13 

D 14 

D21 

D22 

D23 

D24 

D31 

D32 

D33 

-D34 

D41 

Da2 

D43 

D44 

(2.21) 


as  follows: 


GDOP  —  \J  Du  +  D22  +  -D33  +  -D44  (2.22) 

The  GDOP  relates  User  Equivalent  Range  Error  (UERE)  with  root-sum-square  of 
errors  as: 


0?  +  <  +  <  +  <,  =  GDOP  aUERE  (2.23) 

where  ctuere  is  the  standard  deviation  of  UERE. 
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There  are  other  DOP  parameters  which  are  useful  to  characterize  the  accuracy 
of  various  components  of  the  position/time  solution.  These  are  Position  Dilution 
of  Precision  (PDOP),  Horizontal  Dilution  of  Precision  (HDOP),  Vertical  Dilution  of 
Precision  (VDOP),  and  Time  Dilution  of  Precision  (TDOP).  The  definitions  of  these 
DOP  parameters  are  as  follows: 


PDOP 

— 

a/  Dll  +  Do2 

+  d33 

(2.24) 

HDOP 

= 

\J  Du  +  D22 

VDOP 

= 

V  -D33 

TDOP 

= 

V  D44/c 

and  the  relationships  with  c>uere 

are  as  follows: 

\Ae  + 

+  &1 

=  PDOP 

&UERE 

(2.25) 

\f°l 

+  an 

=  HDOP 

&UERE 

=  VDOP 

&UERE 

® Stu 

=  TDOP 

O'UERE 
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2.5  Genetic  Algorithms 


Genetic  Algorithms  (GAs)  are  adaptive  methods  which  are  used  to  solve  search 
and  optimization  problems.  GAs  are  based  on  the  genetic  processes  of  biological 
organisms.  Natural  populations  evolve  over  many  generations  according  to  the  prin¬ 
ciples  of  natural  selection  and  survival  of  the  fittest.  When  they  are  suitably  encoded, 
GAs  are  able  to  evolve  solutions  to  real  world  problems  by  mimicking  the  biological 
evolutionary  process  [6]. 

GAs  operate  on  a  population  of  potential  solutions  applying  the  principle  of  sur¬ 
vival  of  the  fittest  to  generate  better  approximations  to  a  solution.  At  each  generation, 
a  new  set  of  potential  solutions  is  generated  by  the  process  of  selecting  individuals 
according  to  their  level  of  fitness  in  the  population  and  breeding  them  using  operators 
borrowed  from  natural  genetics.  This  process  leads  to  the  evolution  of  populations 
of  individuals  which  have  better  levels  of  fitness  than  the  individuals  that  they  were 
created  from  [34], 

The  implementation  of  GAs  begins  with  an  initial  set  of  parameter  combina¬ 
tions.  Each  individual  of  the  population  is  evaluated  by  a  fitness  function.  A  selection 
process  is  used  to  select  a  group  of  parents  among  the  individuals  with  the  highest 
fitness  values  in  the  population.  The  resulting  parent  chromosomes  (individuals)  are 
paired  off.  The  string  content  is  swapped  between  the  parents  in  a  pair  during  a 
crossover  process,  according  to  a  specified  crossover  probability.  The  resulting  chro¬ 
mosomes  are  called  children  [5].  The  crossover  operation  is  shown  in  Figure  2.4. 
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Crossover  point 


Figure  2.4:  Crossover  operator  cuts  both  chromosomes  at  a  determined  point  and 
swaps  the  tails. 

A  mutation  operation  is  applied  to  each  child  after  crossover.  It  basically  alters 
each  gene  with  a  small  probability  [6].  An  example  of  mutation  is  shown  in  Figure 
2.5. 
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1 
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1 

1 

Figure  2.5:  Mutation  operator  flips  the  value  of  the  bit  at  determined  point. 


After  selection,  crossover,  and  mutation  operations  applied  to  the  current  gen¬ 
eration. ,  a  whole  new  generation  is  produced.  This  new  generation  inherits  a  higher 
proportion  of  the  characteristics  possessed  by  the  good  members  of  the  previous  gen¬ 
eration.  In  this  way,  good  characteristics  are  spread  throughout  the  population  over 
many  generations.  By  favoring  the  mating  of  the  more  fit  individuals,  the  most 
promising  regions  of  the  search  space  are  explored  to  find  the  optimum  solution.  If 
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the  GA  has  been  well  designed  for  the  nature  of  the  problem,  the  population  con¬ 
verges  to  an  optimal  or  a  close  solution  to  the  problem  [6].  This  process  is  shown  in 
Figure  2.6: 


C“D 


Figure  2.6:  Basic  GA  Flowchart 

2.5.1  Selection.  Selection  involves  choosing  the  individuals  in  the  popula¬ 
tion  that  will  create  children  for  the  next  generation.  The  purpose  of  the  selection 
is  to  emphasize  the  fitter  individuals  in  the  population  in  hopes  that  their  children 
will  in  turn  have  even  better  fitness.  Selection  has  to  be  balanced:  too-strong  selec¬ 
tion  means  that  sub-optimal  highly  fit  individuals  takes  over  the  population,  reducing 
the  diversity  needed  for  further  progress;  too-weak  selection  results  in  slow  conver¬ 
gence.  Numerous  selection  schemes  have  been  proposed  in  the  GA  literature  [22], 
The  roulette  wheel  selection  and  tournament  selection  are  explained  briefly. 

A  selection  mechanism  in  GA  is  simply  a  process  that  selects  better  individuals 
in  the  population  for  the  mating  pool.  The  selection  pressure  is  the  degree  to  which 
better  individuals  are  selected:  the  higher  the  selection  pressure,  the  more  the  better 
individuals  are  selected.  This  selection  pressure  drives  the  GA  to  improve  the  pop¬ 
ulation  fitness  over  succeeding  generations.  The  convergence  rate  of  a  GA  is  largely 


determined  by  the  selection  pressure.  If  the  selection  pressure  is  too  low,  the  conver¬ 
gence  rate  will  be  slow,  which  causes  the  GA  to  unnecessarily  take  longer  to  find  the 
optimal  solution.  If  the  selection  pressure  is  too  high,  there  is  an  increased  probability 
for  the  GA  to  converge  to  an  incorrect  (sub-optimal)  solution  [21]. 

2. 5. 1.1  Roulette  Wheel  Selection.  This  approach  is  the  commonly  used 
type  of  fitness-proportional  selection.  With  the  roulette  wheel  selection  approach,  the 
probability  of  selection  is  proportional  to  an  individual’s  fitness  value.  The  selection 
probability  is  the  rate  of  individual’s  fitness  to  the  sum  of  the  population’s  fitness. 
With  this  approach,  fitter  individuals  have  more  probability  of  being  selected.  How¬ 
ever,  even  the  fittest  may  not  be  selected  for  parenthood  [9]. 

2.5. 1.2  Tournament  Selection.  The  tournament  selection  is  made  by 
taking  s  competitors  from  the  population  and  selecting  the  the  fittest  one,  with  s  being 
the  tournament  size.  Therefore,  the  selected  individuals  have  a  better  average  fitness 
than  the  population.  This  difference  drives  the  GA  to  improve  the  fitness  of  each 
succeeding  generation.  Increased  number  of  competitors  causes  increased  probability 
of  early  convergence  [21]. 

2.5. 1.3  Elitism.  Fitness-proportional  selection  does  not  guarantee 
the  selection  of  fitter  individuals,  including  the  fittest.  Therefore,  the  best  solution 
to  the  problem  discovered  might  be  thrown  away.  This  is  sometimes  advantageous, 
because  that  causes  the  algorithm  to  explore  more  search  space  before  convergence. 
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For  many  applications,  the  search  performance  can  be  greatly  improved  by  keeping 
the  best  individual  between  generations.  Designing  the  algorithm  so  that  the  best,  or 
elite,  will  not  die  is  termed  elitism.  In  elitism,  the  elite  is  kept  undisrupted  between 
generations  while  it  is  the  best  [9]. 

2.5.2  Exploration  and  Exploitation.  In  order  to  find  a  global  minimum,  an 
efficient  optimization  algorithm  must  use  two  techniques:  exploration  to  investigate 
new  and  unknown  areas  in  the  search  space,  and  exploitation  to  make  use  of  knowledge 
found  at  points  previously  visited  to  help  find  better  points.  Since  these  requirements 
are  contradictory,  a  good  search  algorithm  must  find  a  tradeoff  between  exploration 
and  exploitation.  A  good  combination  of  these  strategies  can  be  quite  effective,  but 
it  is  hard  to  find  the  balance  between  them  to  find  an  acceptable  solution  [5]. 

2. 6  Summary 

Brief  overview  of  the  GPS  is  given  in  the  beginning  of  this  chapter.  The  out¬ 
lines  of  position  estimation  and  Dilution  of  Precision  concept  in  GPS  are  discussed. 
The  fundamental  concepts  and  operations  for  GA  are  introduced  and  explained.  In 
Chapter  III,  the  methodology  of  the  thesis  is  explained. 
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III.  Methodology 


3. 1  Overview 

The  purpose  of  this  chapter  is  to  describe  the  development  of  the  Matlab® 
simulation  for  this  study.  In  Section  3.2,  the  structure  of  the  simulation  is  provided 
to  give  a  big  picture  of  the  process.  In  Section  3.3,  the  implementation  of  the  GA  and 
reason  for  each  adopted  approach  are  explained.  In  Section  3.4,  the  formulas  used  for 
orbit  propagation  are  given.  In  Section  3.5,  DOP  calculation  is  explained.  In  Section 
3.6,  the  cost  function  to  evaluate  the  candidate  solutions  is  introduced. 

3.2  The  Structure  of  the  Simulation 

The  simulation  is  developed  in  Matlab®  environment  for  this  study  [1].  The 
simulation  includes  a  Genetic  Algorithm  module,  an  orbit  propagator  module,  and 
a  performance  evaluation  module.  The  Genetic  Algorithm  module  generates  candi¬ 
date  solutions  and  provides  them  to  orbit  propagator  module.  The  orbit  propagator 
calculates  the  satellite  positions  over  time.  The  performance  evaluation  module  ex¬ 
ploits  the  satellite  positions  calculated  by  the  orbit  propagator  module  and  the  given 
receiver  positions  to  find  out  the  cost  of  each  candidate  solution.  The  cost  is  then 
provided  to  Genetic  algorithm  module.  The  Genetic  Algorithm  module  assesses  the 
cost  of  each  candidate  solution  and  generates  new  candidates.  This  process  repeats 
until  the  predetermined  criteria  are  met.  This  process  is  illustrated  in  Figure  3.1. 
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Cost  of  Each  Individual 


Figure  3.1:  Top-Level  Simulation  Structure 

3.3  Genetic  Algorithm  Implementation 


Since  there  is  no  unique  genetic  algorithm  that  works  for  all  types  of  opti¬ 
mization  problems,  a  special  algorithm  has  been  developed  for  this  study.  Since  the 
optimization  problem  of  this  study  is  a  minimization  problem,  using  the  term  cost  is 
more  sensible  than  fitness.  An  initial  population  is  generated  randomly,  and  cost  for 
each  individual  is  calculated.  The  individuals  are  mated  randomly  and  the  individu¬ 
als  selected  by  a  tournament  selection  as  parents  for  crossover  operation  to  generate 
an  intermediate  population.  After  the  crossover  operation,  the  mutation  operation  is 
applied  to  the  intermediate  population  to  create  the  children  population.  The  cost 
for  each  child  is  calculated.  A  second  tournament  selection  is  made  over  the  initial 
(old)  population  and  the  children  generation  completely  to  create  the  new  genera¬ 
tion.  This  process  repeats  until  the  stopping  criteria  for  the  algorithm  are  met.  The 
implementation  of  GA  for  this  study  is  illustrated  in  Figure  3.2. 
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Figure  3.2:  Process  for  Generating  New  Populations  (Also  Referred  to  as  a  New 
Generation) 


3.3.1  Representation. 


3. 3. 1.1  Encoding.  The  GA  developed  for  this  study  relies  on  binary 
encoding,  which  is  the  traditional  GA  approach  [13].  Providing  the  required  quan¬ 
tization  level  for  desired  resolution  of  parameters,  the  binary  encoding  is  capable  of 
creating  a  large  enough  search  space. 

In  order  to  encode  the  parameter  values  in  binary  format,  the  number  of  the 
required  bits  should  be  determined  using  Equation  3.1  [13],  where  Zg  is  the  chromosome 
length  (number  of  bits)  of  the  parameter,  t;min  is  the  minimum  parameter  value,  (mai 
is  the  maximum  parameter  value  and  R  is  the  desired  resolution  (number  of  digits 
after  the  decimal  point). 
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<  (f 


Uin)10fi  <  2‘f  -  1 


(3.1) 


Calculating  the  required  chromosome  length  for  the  parameter  value  allows  to 
represent  the  parameter  value  in  the  Matlab  code  for  GA  operations.  Then  the  binary 
coded  parameter  value  (Cbm)  can  be  converted  to  its  decimal  value  ( Cdec )•  The  decimal 
value  (Cdec)  should  be  converted  to  its  adjusted  or  real  value  (Creai)  by  using  Equation 
3.2  for  calculations. 


j-  _  £  Craax  Crain  £ 

£,real  s  min  i  ^  S^e 


(3.2) 


3.3. 1.2  Population  Size.  It  is  very  important  to  create  an  initial  pop¬ 
ulation  which  provides  a  good  trade-off  between  efficiency  and  effectiveness.  When 
the  population  is  small,  the  algorithm  may  not  effectively  explore  the  search  space. 
On  the  other  hand,  when  the  population  is  large,  the  efficiency  of  the  algorithm  is 
significantly  reduced  to  find  an  optimal  solution  in  a  reasonable  time  of  computa¬ 
tion.  Therefore,  there  should  be  an  optimal  population  size  for  a  given  individual  bit 
length.  At  least,  a  minimum  population  size  should  be  determined  which  enables  the 
algorithm  to  visit  every  point  in  the  search  space.  That  minimum  population  size  can 
be  calculated  using  Equation  3.3  [25]. 


n  =  (1  -  (1/2 )N™n-iy 


(3.3) 
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Where  II  is  the  probability  of  points  in  the  search  space  to  be  explored  and  Nmin 
is  the  minimum  population  size.  Using  an  exponential  function  approximation  the 
equation  can  be  written  as  follows: 


(3.4) 


]T[  r^,  g  2^ min  ■*- 


Then  minimum  population  size  can  be  calculated  using  the  following  equation: 


(3.5) 


In  this  study,  the  population  size  is  determined  so  that  it  is  not  less  than  the 
minimum  population  size  and  it  is  not  too  large  (which  decreases  the  efficiency  of  the 
algorithm). 

3.3.2  Selection.  In  Holland’s  original  genetic  algorithm,  parents  are  replaced 
by  their  children  just  after  they  give  birth  which  is  called  generational  replacement. 
Because  children  may  be  worse  than  their  parents,  some  fitter  chromosomes  will  be 
lost  from  the  evolutionary  process  with  the  strategy  of  replacing  each  parent  with  his 
children  directly.  Several  strategies  have  been  developed  to  overcome  this  problem 
[13].  In  this  thesis  research,  most  of  the  generated  chromosomes  have  poor  fitness 
values;  in  fact,  some  chromosomes  have  fitness  values  of  not-a-number  (NaN)  if  there 
is  no  solution.  Therefore,  there  is  no  luxury  of  losing  good  chromosomes.  By  applying 
Holland’s  original  genetic  algorithm  and  his  generational  replacement,  the  population 
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of  the  problem  in  this  study  was  not  able  to  converge  properly.  In  order  to  have  an 
evolving  population,  enlarged  sampling  space  approach  has  been  adopted  in  this  study. 
When  selection  is  performed  on  enlarged  sampling  space,  both  parents  and  children 
have  the  same  chance  of  competing  for  survival.  This  approach  evidently  has  an 
advantage,  in  that  genetic  algorithm  performance  can  be  improved  by  increasing  the 
crossover  and  mutation  rates.  There  is  no  need  to  worry  that  the  high  rate  introduces 
too  much  random  perturbation,  if  selection  is  performed  on  the  enlarged  sampling 
space  [13]. 

There  are  many  selection  schemes  for  GAs,  each  with  different  characteristics. 
An  ideal  selection  scheme  would  be  simple  to  implement.  It  also  should  be  able  to 
adjust  its  selection  pressure  to  tune  its  performance  for  different  domains.  Tourna¬ 
ment  selection  increasingly  being  used  as  a  GA  selection  scheme,  because  it  satisfies 
above  mentioned  criteria. 

Tournament  selection  was  used  in  the  algorithm  due  to  its  resistance  to  early 
convergence.  The  tournament  selection  approach,  which  closely  mimics  mating  com¬ 
petition  in  nature,  is  to  randomly  pick  a  small  subset  of  chromosomes  (two  or  three) 
from  the  mating  pool,  and  the  chromosome  with  the  lowest  cost  in  this  subset  becomes 
a  parent.  The  tournament  repeats  for  every  parent  needed.  Tournament  selection 
works  best  for  larger  population  sizes,  because  sorting  becomes  time-consuming  for 
large  populations  [14]. 
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The  tournament  selection  with  two  chromosomes  in  this  study  ensures  that  the 
best  individual  does  not  die  and  the  selection  pressure  is  lowest  possible  to  avoid 
excessive  elitism. 

Tournament  selection  is  made  twice  in  a  generation  in  this  study.  The  first  one 
is  made  to  select  parents  for  producing  intermediate  population;  A  number  of  couples 
are  selected  out  of  N  number  of  parents,  where  A  is  the  number  of  children  population 
size  and  N  is  the  population  size.  The  second  tournament  selection  is  accomplished  to 
produce  the  next  generation  from  the  current  population;  N  number  of  new  parents 
are  selected  out  of  N  +  A  number  of  parents  and  children  of  current  population.  At 
each  time  two  random  individuals  are  selected  and  their  fitness  values  are  compared 
and  better  one  is  selected. 

3.3.3  Crossover.  After  having  selected  the  mates  for  crossover,  the  ques¬ 
tion  arises:  How  many  crossovers  should  be  made?  There  are  two  issues  relating  to 
crossover:  bias  and  disruption.  The  investigation  of  biasing  effect  of  crossover  showed 
that  there  are  two  sources  of  bias  exist  to  be  exploited  in  a  genetic  algorithm:  posi¬ 
tional  bias,  and  distributional  bias.  One-point  crossover  has  considerable  positional 
bias,  which  causes  early  convergence  to  sub-optimal  solution.  On  the  other  hand, 
one-point  crossover  has  no  distributional  bias,  because  the  crossover  point  is  chosen 
randomly.  However,  lack  of  this  bias  is  not  always  good  since  it  limits  the  exchange 
of  information  between  the  parents.  A  completely  random  number  crossover,  which 
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is  called  uniform  crossover ,  is  capable  of  significantly  reducing  any  biases.  In  uniform 
crossover,  each  crossover  operation  is  made  in  random  numbers  [25]. 

Uniform  crossover  is  also  capable  of  improving  the  performance  of  algorithm 
as  explained  by  Spears  and  De  Jong  [29].  However,  uniform  crossover  causes  disrup¬ 
tion.  Because  the  mates  are  cut  multiple  times,  the  good  features  of  the  parents  are 
disrupted  so  that  the  children  are  less  likely  to  have  good  fitness  values.  Disruption 
slows  down  evolution  in  the  beginning  of  the  algorithm,  but  it  may  be  helpful  as  the 
algorithm  proceeds.  When  the  population  becomes  quite  homogenous,  the  ability 
of  crossover  to  produce  new  individuals  is  reduced.  This  ability  is  called  crossover 
productivity.  Uniform  crossover  is  capable  of  having  more  productivity  as  evolution 
proceeds.  Uniform  crossover  is  used  in  this  study  because  of  these  factors. 

In  order  to  implement  uniform  crossover,  a  randomly  generated  vector,  which 
is  called  mask,  at  the  same  length  of  the  individuals  is  used.  In  the  mask  vector,  a 
switch  from  0  to  1  or  from  1  to  0  represents  a  crossover  point  at  respective  bit  of 
mates.  This  process  can  be  formulated  as  follows: 

01  =  m  0  (5i  ©  m  <8)  f32  (3.6) 

02  —  Jn  0  /3i  ©  m  (g) 

where  (3\  and  f3 2  are  parents,  m  is  mask  vector,  m  is  complement  of  mask  vector, 
0X  and  02  are  children  and  ©,  0  denote  component-wise  addition  and  multiplication 
respectively  [25].  This  process  is  illustrated  in  Figure  3.3. 
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Figure  3.3:  In  uniform  crossover,  a  switch  from  0  to  1  or  from  1  to  0  represents  a 
crossover  point  at  respective  bit  of  mates. 

3.3.4  Mutation.  Mutation  operation  operation  is  accomplished  by  determin¬ 
ing  the  number  of  bits  to  be  mutated  in  the  population.  This  number  is  calculated  by 
multiplying  the  probability  of  mutation  by  the  population  size  and  the  chromosome 
length  as  shown  in  Equation  3.7  [14]. 


nm  =  PmNl  (3.7) 

where  nm  is  the  number  of  bits  to  be  mutated,  Pm  is  mutation  probability  and  l  is 
total  chromosome  length.  Then  nm  bits  are  randomly  flipped  in  the  population.  This 
process  is  illustrated  in  Figure  3.4. 

3.4  Orbit  Propagation 

The  orbit  propagation  module  requires  six  orbital  elements  to  calculate  the 
position  of  satellites.  These  elements  and  their  definitions  are  as  follows: 
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Figure  3.4:  In  distributed  mutation  operation,  number  of  bits  to  be  mutated  (nm) 
are  calculated  and  nm  bits  are  randomly  flipped  in  the  generation. 


a  Semi-major  axis,  defines  the  size  and  period  of  the  orbit 


e  Eccentricity,  defines  the  shape  of  the  orbit 


i  Inclination,  defines  the  angle  between  orbit  plane  and  the  equator 

u  Argument  of  perigee,  defines  the  angle  between  the  ascending  node  and  the  orbit’s 
point  of  closest  approach  to  the  earth 


v  True  anomaly,  defines  the  angle  between  perigee  and  the  satellite  position  vector 

Right  Ascension  of  Ascending  Node  (RAAN),  defines  the  angle  between  vernal 
equinox  and  the  point  where  the  orbit  crosses  the  equatorial  plane  from  south 
to  north 


In  Figure  3.5  [4],  these  elements  are  illustrated. 
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perigee  direction 


(b) 

Figure  3.5:  (a)  Three  dimensional  view  of  the  orbital  elements, 

(b)  Orbital  elements  defined  within  the  orbital  plane  [4], 
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In  this  study,  using  Longitude  of  Ascending  Node  (LAN)  (which  is  defined 
relative  to  the  earth)  is  more  sensible  than  RAAN  (which  is  defined  relative  to  the 
inertial  space),  because  the  satellites  of  the  objective  constellation  should  be  visible 
over  a  specific  region  of  the  earth.  The  input  for  the  orbit  propagation  is  Longitude  of 
Ascending  Node,  but  the  corrected  RAAN  is  used  for  calculation  of  the  actual  orbit. 
This  calculation  ensures  that  the  ground  track  of  the  satellite  will  be  same  for  different 
true  anomaly  values  where  the  other  parameters  are  constant.  Therefore,  it  will  be 
possible  to  make  an  optimization  of  true  anomaly.  In  Figure  3.6,  the  ground  tracks 
of  Satl  (red)  and  Sat2  (blue)  are  shown,  which  have  the  same  orbital  parameters 
except  the  true  anomaly.  The  true  anomaly  values  are  180°  and  210°  respectively. 
Because  the  RAAN  is  not  corrected  to  reflect  the  LAN,  the  ground  track  of  Sat2  drifts 
approximately  30°  towards  east.  The  angular  difference  of  two  ground  tracks  on  the 
equatorial  plane  is  the  required  correction  to  reflect  the  same  LAN  value.  Therefore, 
the  corrected  RAAN  is  calculated  as  follows: 

^  =  ^LAN  —  Qshift  (3.8) 

where  FIlan  is  Longitude  of  Ascending  Node  and  shift  is  the  required  correction  for 
RAAN. 

Argument  of  perigee  has  a  similar  effect  on  the  ground  track  too.  Therefore, 
the  effect  of  both  argument  of  perigee  and  true  anomaly  should  be  calculated.  The 
value  of  total  effect  can  be  calculated  as  follows: 
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Figure  3.6:  The  Effect  of  the  Uncorrected  RAAN.  If  the  RAAN  is  not  corrected  to 
represent  LAN,  the  ground  track  drifts  when  the  true  anomaly  is  changed. 


Q 


Mu  +  Mv 


shift 


n 


+  &GMST 


(3.9) 


where  Mu  is  the  mean  anomaly  equivalent  of  argument  of  perigee,  Mv  is  the  mean 
anomaly  equivalent  of  true  anomaly,  n  is  mean  motion,  Vte  is  the  earth’s  rotation 
rate,  and  Ogmst  is  Greenwich  mean  sidereal  time.  9gmst  can  be  assumed  zero  for 
simplification  and  mean  motion  is  defined  as: 


(3.10) 

where  /i  is  the  earth’s  gravitational  constant.  Before  calculating  the  mean  anomaly 
equivalent  of  true  anomaly  and  argument  of  perigee,  their  eccentric  anomaly  equiva¬ 
lents  should  be  calculated. 
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tan(-^ 


(3.11) 


.  / 1  —  e  u 

)  =  y—et&n{2] 

tanT'j  =  \/yT7  tan(|)  (3-12) 

where  is  the  eccentric  anomaly  equivalent  of  argument  of  perigee  and  ly,  is  eccen¬ 
tric  anomaly  equivalent  of  true  anomaly.  Then  the  equivalent  mean  anomaly  values 
are  calculated  as  follows: 


Mu  =  Eu-  e  sin (EJ) 


(3.13) 


Mv  =  Eu  —  e  sin (Ev) 


(3.14) 


The  position  of  the  satellites  are  calculated  a  for  period  of  simulation  time  (tsim) 
at  simulation  steps  {tstep). 


Mstep  =  nt  +  Mv 


(3.15) 


Estep  -Mstep  “I-  ^ 


(3.16) 
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(3.17) 


tan( 


^ step 


) 


tan(^f) 


Now,  the  position  and  velocity  vectors  can  be  calculated  using  two-body  equa¬ 
tions  as  follows  [33]: 


p  =  a  (1  —  e2) 


(3.18) 


YPQW 


P  COS [Ustep] 

l+ecos(i/steP) 

p  sin (vstep) 

l+ecos(ystep) 

0 


(3.19) 


VPQW 


\J~^p  sin(I/step) 

(e  +  cos(z/step)) 

0 


(3.20) 


where  p  is  orbit  semiparameter,  rpQ\y  is  position  vector  in  perifocal  coordinate  system 
and  v pqw  is  velocity  vector  in  perifocal  coordinate  system.  It  is  required  to  convert 
those  vectors  from  perifocal  coordinate  system  to  ECI  frame.  The  following  rotation 
matrices  should  be  calculated  to  do  this  conversion: 
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R 


cos( — uj')  sin(— oj)  0 


—  sin(— a;)  cos(— u;)  0 


(3.21) 


R; 


R 


a  — 


0  cos  (—i)  sin(— i) 

0  —  sin(— i)  cos( — f ) 

cos(— il)  sin(— f2)  0 

—  sin(— Q)  cos(— Q)  0 


0 


0 


(3.22) 


(3.23) 


Using  the  above  matrices,  the  position  and  velocity  vectors  can  be  converted  as  follows: 


r ec i  —  Rn  R)  Ru  rpQw 


(3.24) 


v ec i  —  Rn  R/  Rl,  vpqw 


(3.25) 


Coordinate  frame  conversion  is  not  complete  yet.  Conversion  from  ECI  frame 
to  ECEF  frame  is  also  required.  In  order  to  convert  ECI  frame  to  ECEF  frame,  how 
far  the  earth  rotated  in  tstep  seconds  should  be  found. 
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(3.26) 


9  Ue  tstep 

where  6  is  the  earth’s  rotation  angle.  The  following  rotation  matrix  can  be  calculated 
using  8\ 


T  EC  ItoECEF 


cos (9)  sin(6?)  0 

—  sin(6l)  008(6*)  0 
0  0  1 


(3.27) 


where  T ECitoECEF  is  rotation  matrix  from  ECI  to  ECEF  frame.  Using  this  rotation 
matrix,  the  position  and  velocity  vectors  can  be  converted  from  ECI  to  ECEF  frame 
as  follows: 


r ECEF  —  T ECitoECEF  r ECI  (3.28) 

v ECEF  =  T ECitoECEF  V ECI  (3.29) 

The  position  vector  of  a  satellite  in  ECEF  frame  ( tecef )  calculated  in  Equation 
3.28  is  needed  to  calculate  DOP  value.  In  the  simulation  two-body  equations  of 
motion  are  used  to  optimize  the  constellation.  The  accuracy  of  two-body  equations 
are  enough  for  short  term  simulations.  However,  two-body  equations  are  not  realistic 
for  long  term  simulations,  where  higher  order  terms  have  a  larger  impact.  Therefore, 
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J2  perturbation  is  included  in  computation  for  long  term  analysis  of  the  solutions. 
In  order  to  include  J2  perturbation  effect,  the  initial  condition  of  satellites  should 
be  calculated  using  Equation  3.11  through  Equation  3.25.  The  calculation  of  initial 
position  and  velocity  vectors  enables  to  form  the  current  6x1  state  matrix.  One  can 
calculate  the  next  state  including  the  J2  perturbation  as  follows: 


" 

0x3 

= 

03x3 

hx3 

0x3 

+ 

03x3 

03x3 

0x3 

Vlx3 

633,3 

03x3 

Vlx3 

a3x3 

03x3 

Vlx3 

t+ 1 

t 

- 

(3.30) 


where  033,3  is  zero  matrix  and  1 3^3  is  identity  matrix.  633,3  and  a33,3  are  defined  as 
follows: 


(3/0  0 

0  aj  0 

0  0  CLx 


(3.31) 


(3.32) 


<3/,  aj  and  a,K  are  the  accelerations  due  to  J2  perturbation  in  the  respective  axis  in 
ECI  frame.  J2  acceleration  components  are  defined  as  follows  [33]: 
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ai  = 


2 


(3.33) 


(1 


5|M|2 


2 


(3.34) 


2 


(3.35) 


where  R®  is  earth’s  radius.  After  the  satellite  positions  are  calculated  in  ECI  frame, 
they  are  converted  to  ECEF  frame  using  the  equations  3.26,  3.27  and  3.28. 

3.5  DOP  Calculation 

The  DOP  should  be  calculated  for  different  receiver  positions  in  the  desired 
region.  The  receiver  position  coordinates,  which  are  defined  in  the  World  Geodetic 
System  1984  (WGS84)  reference  system,  should  be  converted  to  ECEF  frame.  In 
this  study,  this  conversion  is  accomplished  using  the  Ua2ecef  function  of  Aerospace 
toolbox  in  MATLAB.  After  having  the  receiver  positions  converted  to  ECEF  frame, 
one  can  determine  the  elevation  angle  of  the  Line-of-Sight  (LOS)  vectors  from  receiver 
to  satellites  using  the  elevation  function  of  Mapping  toolbox  in  MATLAB.  Then 
elevation  angles  (El)  for  each  satellite  should  be  more  than  mask  angle  ( Elmask ): 


(3.36) 
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The  satellites  which  are  below  mask  angle  are  assumed  as  iron-visible.  If  there  are 
fewer  than  four  visible  satellite,  DOP  is  undefined  and  valued  as  NaN.  If  there  are  at 
least  four  satellites,  the  unit  LOS  vectors  are  calculated  as  follows: 


uk  =  rk  -  ru  (3.37) 

where  U*,  is  the  receiver  to  satellite  LOS  unit  vector,  rk  is  position  vector  of  kth  satellite 
in  ECEF  frame  (calculated  as  yEcef  in  the  previous  section)  and  ru  is  receiver  position 
vector  in  ECEF  frame.  Then  the  H  is  formed  as  follows: 


H 


Lp  —  1 
U2  -1 


(3.38) 


The  DOP  matrix  is  calculated  using  Equation  2.21.  Because  GDOP  represents 
the  overall  satellite  geometry,  GDOP  value  is  calculated  to  analyze  the  navigation  per¬ 
formance  of  the  constellation  using  Equation  2.22.  The  system  is  considered  available 
if  the  GDOP  value  is  less  than  a  threshold  value  ( GDOPt ).  Therefore,  the  system 
must  meet  the  following  criterion: 


GDOP  <  GDOPt 


(3.39) 
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After  calculating  GDOP  value  and  availability  for  all  receivers  and  time  steps, 
the  proposed  solution  can  be  evaluated.  In  this  study,  performance  evaluation  of  the 
constellation  is  made  considering  accuracy  and  availability. 

3. 6  Cost  Calculation 

The  GA  developed  for  this  study  requires  a  single  cost  value  for  each  candidate 
solution.  The  cost  function  is  developed  to  evaluate  the  candidate  solutions  based  on 
the  following  criteria: 

1.  The  availability  should  be  100%,  if  possible. 

2.  Average  GDOP  should  be  the  lowest  possible. 

3.  The  system  should  perform  better  over  the  region  of  interest. 

The  cost  function  developed  to  meet  the  above  requirements  is  as  follows: 

Tlr  Tl'p 

Cf  =  Wqdop  E  Wi  GDOP i  +  Wa  Wi  (1  -  on)  (3.40) 

2=1  2=1 

where  Cf  is  cost  function,  Wqdop  is  GDOP  weight,  nr  is  number  of  the  receiver 
points,  Wi  is  weight  of  ith  receiver  point,  GDOPi  is  average  GDOP  value  at  the  ith 
receiver,  Wa  is  availability  weight,  and  a*  is  availability  at  the  ith  receiver  (expressed 
as  a  decimal  notability,  where  1  is  always  available  and  0  is  never  available).  NaN 
GDOP  values  are  disregarded  for  average  GDOP  calculation. 
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For  most  applications,  which  are  using  satellite  based  navigation,  it  is  desired 
to  have  the  system  at  all  times  with  the  most  accuracy.  Therefore,  the  navigation 
system  should  provide  a  good  accuracy  with  continuous  service.  The  cost  function 
was  developed  to  meet  this  requirement. 

3. 7  Summary 

In  this  chapter,  the  structure  of  the  simulation  developed  for  this  study  is  cov¬ 
ered.  The  adopted  approach  for  each  operator  of  GA  is  explained  in  detail.  The 
equations  of  the  orbit  propagator  and  DOP  calculation  are  given.  Finally,  the  cost 
evaluation  function  is  formed  to  distinguish  the  better  constellations  among  the  pop¬ 
ulation.  In  Chapter  IV,  the  obtained  data  and  analysis  are  given. 
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IV.  Data  and  Analysis 


4-1  Overview 

In  this  chapter,  the  obtained  data  are  given  and  analyzed.  I11  Section  4.2,  the 
settings  of  the  developed  simulation  are  given.  In  Section  4.3,  the  nature  of  the  search 
space  of  this  problem  is  analyzed.  I11  Section  4.4,  the  designed  GA  is  tested  to  see  if  it 
is  effective  for  generating  solutions  for  the  problem  of  this  thesis  effort.  In  Section  4.5, 
the  obtained  data  for  the  cases  are  analyzed  in  groups,  which  have  similar  constraints. 
Finally,  in  Section  4.6,  the  J2  perturbation  effect  on  the  designed  constellations  with 
different  features  are  analyzed. 

4-2  Simulation  Settings 

All  but  one  constellation  proposed  in  this  thesis  research  are  comprised  of  five 
satellites.  Several  cases  are  designed  with  different  constraints.  All  cases  were  run 
for  three  different  random  seeds  to  mitigate  the  bias  caused  by  Matlab®’s  random 
number  generator.  One  of  the  cases  was  designed  to  have  four  satellites. 

Based  on  the  constraints  determined  for  the  cases  and  the  resolution  of  the 
variables,  total  bit  length  for  each  case  were  calculated  according  to  using  Equation 
3.5.  The  resolution  of  the  variables  and  the  required  bit  length  calculation  for  Case  1 
are  given  in  Table  4.1. 

The  sum  of  the  bit  length  of  each  variable  for  the  first  case  equals  47  for  one 
satellite  in  Case  1,  which  is  the  unconstrained  case.  Then  the  total  bit  length  (7)  is 
calculated  as  (47x5  =  235).  The  next  step  is  to  find  the  minimum  population  size. 
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Table  4.1:  The  required  bit  length  calculation  for  the  first  case 


e 

i 

^LAN 

C 0 

V 

R 

2 

1 

1 

0 

0 

£ min  •  £,max 

0  :  0.85 

0  :  90 

0  :  360 

0  :  360 

0  :  360 

k 

7 

10 

13 

9 

9 

In  order  to  ensure  almost  a  100%  probability  of  the  algorithm  to  visit  every  point 
of  the  search  space  with  crossover  operation  and  without  mutation,  the  minimum 
population  size  for  l  =  235  is  calculated  as  55.  In  the  simulation  for  every  case,  the 
population  size  is  used  as  400  which  is  more  than  the  minimum  required  population 
size  and  provides  a  better  algorithm  performance.  The  number  of  the  generations 
for  the  cases  was  set  to  300,  600,  and  900,  depending  on  the  convergence  rate  of  the 
cases. 

Case  1  has  (2235  =  5.52x10'°)  candidate  solutions  and  is  run  with  a  parent 
population  number  of  400  and  children  number  of  400  for  600  generations.  Calculating 
the  cost  of  every  candidate  solution  would  be  impossible,  where  the  GA  calculates 
the  cost  of  only  a  maximum  of  (800x600  =  480000)  candidate  solutions  in  an  effort 
to  find  an  optimal  solution  in  a  reasonable  computation  time.  Obviously,  the  GA 
evaluates  only  very  small  amount  of  the  possible  candidate  solutions.  That  is  one  of 
the  prominent  features  of  the  GA.  Because  it  is  not  possible  to  evaluate  every  possible 
solutions  using  the  GA  with  above  mentioned  parameters,  it  is  not  guaranteed  that 
the  optimal  solution  found  by  the  GA  is  the  best  one.  However,  the  generated  solution 
is  an  optimal  and  reasonable  solution. 
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The  parameters  and  the  constants  used  in  the  simulation  are  given  in  Table  4.2. 


The  constraints  for  each  case  are  given  in  Table  4.3. 


Table  4.2:  The  constants  and  parameters  of  the  simulation. 


Symbol 

Parameter 

Value 

a 

Semi-major  axis 

42164169.6  m 

n 

Earth’s  gravitational  constant 

398600.5  km3/s 2 

P® 

Earth’s  radius 

6378137  m 

tie 

Earth’s  rotation  rate 

7.2921151  xl0~5  rad/s 

J2 

J2  perturbation  coefficient 

0.0010826269 

tsim 

Simulation  time 

259000  s 

t  step 

Simulation  time  step 

900  s 

Elmask 

Elevation  mask  angle 

10° 

GDOPt 

GDOP  threshold  value 

10 

Wgdop 

GDOP  weight 

1 

wa 

Availability  weight 

300 

N 

Parent  Population  size 

400 

X 

Children  population  size 

400 

9 

Number  of  generations 

300,  600,  and  900 

Pc 

Crossover  Probability 

1 

Pm 

Mutation  Probability 

l/l 

55 
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The  simulation  was  run  for  20  cases-18  using  the  GA,  one  random  run  without 


the  GA  and  a  benchmark  run.  The  description  of  the  cases  are  given  in  Table  4.4. 


Table  4.4:  The  Description  of  the  Cases. 

Case  Description 

1  Unconstrained 

2  Four  satellite 

3  2  GEO  DLav(0,70),  3  GSO 

4  2  GEO  (350: 10, 60:80),  3  GSO 

5  2  GEO,  3  GSO 

6  2  GEO  0Lav(330:110),  3  GSO 

7  2  GEO  0Lav(330:110),  3  GSO  *(63.4) 

8  2  GEO  0Lav(330:110),  3  GSO  *(63.4)  u( 270) 

9  2  GEO  0Lav(330:110),  3  GSO  e(0.3:0.7)  *(63.4)  cu(270) 

10  2  GEO  0Lav(330:110),  3  GSO  e(0.3:0.7)  *(63.4) 

11  2  GEO  0^(330:110),  3  GSO  *(63.4)  OiA7V(0:180)  w(270) 

12  2  GEO  Ql/uv(330:110),  3  GSO  *(63.4)  0LA7V(0:180)  w(180:360) 

13  2  GEO  0^(330:110),  3  GSO  e(0.3:0.7)  *(63.4)  QLAN(0:180 )  cu(180:360) 

14  2  GEO,  3  GSO  e(0.3:0.7)  *(63.4)  0LA7V(0:180)  a;(180:360) 

15  2  GEO,  3  GSO  e(0.3:0.7)  *(63.4) 

16  5  GSO  cu(180:360) 

17  5  GSO  *(0:75)  0LAJV(180:360)  o;(180:360) 

18  Unconstrained,  random 

19  2  GEO  OLAiV(0,70),  3  GSO  e(0.5)  *(63.4)  OLAtv(60,80,100)  cu(270)  1/(0:360,0:360,180) 

20  Benchmark  (Same  Constraints  with  Case  19) 


The  constellations  were  optimized  to  have  the  best  performance  within  the  spec¬ 
ified  region  (shown  in  Figure  1.1)  by  applying  the  geographically  distributed  receiver 
weights.  The  coordinates  and  the  weights  of  thirteen  receiver  points  are  given  in  Table 
4.5.  They  are  also  shown  on  the  map  in  Figure  4.1. 
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Tabic  4.5:  The  Coordinates  and  Weights  for  Each  Receiver  Site 


Receiver 

Latitude  (°) 

Longitude  (°) 

Weight 

1 

40  00 

35  00 

5 

2 

43  00 

27  00 

5 

3 

43  00 

43  00 

5 

4 

36  00 

27  00 

5 

5 

36  00 

43  00 

5 

6 

45  00 

20  00 

4 

7 

45  00 

50  00 

4 

8 

35  00 

50  00 

4 

9 

35  00 

20  00 

4 

10 

30  00 

30  00 

3 

11 

30  00 

40  00 

3 

12 

50  00 

30  00 

2 

13 

50  00 

40  00 

2 

Figure  4.1:  Receiver  Weights  and  Geographical  Distributions 
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4-3  The  Nature  of  the  Search  Space 


The  level  of  difficulty  to  find  a  solution  for  this  problem  can  be  inferred  by 
analyzing  the  nature  of  the  search  space.  Because  there  are  25  variables  in  this 
problem,  it  is  hard  to  expose  all  of  the  relationships  between  the  variables.  However, 
a  limited  analysis  can  be  made  by  keeping  all  variables  constant  except  one  variable 
at  a  time. 

The  effect  of  each  variable  on  the  cost  was  calculated  while  keeping  the  other 
variables  constant.  It  is  crucial  to  state  that  this  analysis  does  not  reflect  the  ex¬ 
act  nature  of  the  search  space.  Moreover,  the  results  of  this  analysis  is  specific  to 
the  selected  constellation  (Case  1  Seed  4).  Therefore,  for  an  arbitrary  constellation, 
the  relationships  between  the  variables  might  be  different.  However,  this  analysis  is 
enough  to  give  a  sense  of  the  roughness  of  the  search  space.  The  ground  track  of 
the  constellation  selected  for  this  analysis  is  given  in  Figure  4.7.  Parameters  of  the 
constellation  used  for  this  analysis  is  given  in  Table  4.6. 

The  effect  of  each  variable  on  the  cost  for  the  selected  constellation  is  shown 
in  Figure  4.2.  The  effect  of  the  eccentricity  and  the  inclination  is  not  very  rough. 
However,  the  effect  of  the  longitude  of  ascending  node,  the  argument  of  perigee,  and 
the  true  anomaly  are  very  rough.  There  are  lots  of  peaks  and  dips  in  the  plots  related 
to  those  variables.  This  roughness  of  the  search  space  shows  that  there  are  numerous 
local  minimums,  which  makes  it  hard  to  search  for  global  minimum.  This  situation 
is  more  serious  when  the  variables  are  unconstrained.  Overall,  this  analysis  shows 
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that  the  search  space  of  the  problem  is  very  rough,  it  is  challenging  to  find  a  global 


minimum. 


Table  4.6:  The  Parameters  of  the  Constellation  Used  for  the  Search 
Space  Analysis 


Satellite 

1 

2 

3 

4 

5 

Eccentricity 

0.2075 

0.85 

0.8433 

0.0335 

0.0803 

Inclination 

23.0499 

87.8006 

85.8651 

5.5425 

19.8827 

LAN 

52.4835 

150.6813 

84.6593 

79.8242 

35.7802 

Argument  of  Perigee 

286.7319 

292.3679 

247.9843 

348.0235 

331.1155 

True  Anomaly 

156.3992 

11.9765 

208.5323 

138.0822 

270.5284 

4-4  Validation  of  the  Simulation  Performance 

Before  analyzing  the  solutions  generated  by  the  simulation,  it  should  be  vali¬ 
dated  that  the  GA  is  designed  well  for  the  problem.  Two  tests  were  made  for  this 
purpose.  The  first  one  is  to  test  the  GA  in  a  case  where  the  best  (or  benchmark) 
solution  is  already  known.  The  second  one  is  to  test  if  the  GA  performs  in  a  random 
manner  or  it  evolves. 

4-4-1  Compare  Generated  Solutions  with  a  Benchmark.  Does  the  developed 
simulation  work  well  in  terms  of  GA?  This  question  can  be  answered  by  determining 
the  best  solution  after  evaluating  all  possible  solutions  and  comparing  it  with  the 
solution  generated  by  the  GA. 

As  discussed  before,  it  is  impossible  to  evaluate  every  possible  solution  in  a 
reasonable  amount  of  time.  However,  it  is  possible  to  limit  the  search  space  so  that 
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Cost  Cost 


(a)  (b) 


(c) 


(d) 


(e) 


Figure  4.2:  (a,b,c,d,e,  and  f)  The  effect  of  each  variable  on  the  cost  for  a  selected 

constellation  (e,  i,  Qlan,  oj,  and  u  respectively). 
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every  possible  solution  in  this  limited  search  space  can  be  evaluated  in  a  reasonable 
amount  of  time.  The  benchmark  case,  with  a  limited  search  space,  is  designed  as 
two  GEO  satellites  at  0°  and  70°  on  the  equator  and  three  GSO  satellites  with  e(0.5) 
i(63.4)  0LAiv(60,80,100)  u;(270)  z/(0:360, 0:360, 180).  Only  true  anomaly  vales  (is)  of 
the  third  and  fourth  satellites  are  not  fixed.  Because  true  anomaly  (is)  is  represented 
by  9  bits  (512  quantization  levels  from  0  to  360)  there  are  (512x512=262144)  possible 
solutions  to  be  evaluated.  This  limited  search  space  made  it  possible  to  evaluate  every 
solution  and  find  the  best  one.  Another  case,  Case  19,  was  designed  with  the  same 
search  space  for  the  GA  run.  Case  19  is  run  with  a  parent  population  number  of  50  and 
children  population  number  of  50  for  50  generations.  So,  a  maximum  of  (51x50=2550) 
points  are  visited  at  each  run  for  five  different  seed  numbers.  The  ratio  of  the  visited 
points  to  the  number  of  all  points  in  the  search  space  is  (2550/262144  ss  0.0097). 
It  is  expected  that  the  GA  generates  very  similar  solutions  to  the  solution  found  in 
the  benchmark  case,  which  is  the  best  one.  The  solutions  generated  by  the  GA  and 
the  benchmark  runs  are  given  in  Table  4.7.  The  ground  tracks  of  the  solutions  are 
depicted  in  Figure  4.3. 

The  GA  has  generated  very  similar  solutions  to  the  best  solution.  The  results 
of  this  test  showed  that  the  designed  GA  is  evolving  and  capable  of  generating  good 
solutions.  Therefore,  it  can  be  stated  that  the  GA  has  accomplished  the  expectations. 

4-4 Random  vs  GA.  The  GA  heavily  depends  heavily  upon  randomness 
for  its  results.  Because  the  solution  is  generated  after  evaluation  of  multiple  possible 
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Average  GDOP  :  13.53  Availability  :  72.65%  Cost :  4862.05 


Longitude 

(d) 

Average  GDOP  :  10.56  Availability  :  71.50%  Cost :  4885.03 


Figure  4.3:  (a)  Ground  track  of  the  optimal  solution  for  benchmark  case. 

(b,c,d,e,  and  f)  Ground  tracks  of  the  generated  solutions  for  Case  19  (Random  Seed 
1,2, 3, 4  and  5  respectively).  Constraints:  (2  GEO  L2lav(0,70),  3  GSO  e(0.5)  i(63.4) 
OLatv(60,80,100)  u;(270)  ^(0:360,0:360,180) 
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Table  4.7:  Comparison  of  the  Benchmark  and  Case  19  Runs 


Case 

Average  GDOP 

Availability(%) 

Cost 

A  Cost  (%) 

Benchmark 

12.09 

73.02 

4740.14 

0.00 

Case  19  Seed  1 

10.40 

71.50 

4881.34 

+2.98 

Case  19  Seed  2 

12.97 

73.02 

4773.99 

+0.71 

Case  19  Seed  3 

13.53 

72.65 

4862.05 

+2.57 

Case  19  Seed  4 

10.40 

71.50 

4881.34 

+2.98 

Case  19  Seed  5 

10.56 

71.50 

4885.05 

+3.06 

solutions,  one  might  question  whether  the  GA  performs  any  better  than  a  random 
search.  A  test  was  conducted  to  see  if  it  is  possible  to  find  similar  solutions  to  the 
GA  solutions  by  searching  the  search  space  randomly.  If  it  is  possible,  then  there  is 
no  need  to  use  the  GA.  This  issue  can  be  resolved  by  comparing  the  GA  performance 
and  the  random  search  performance  for  the  same  search  space.  Case  18  was  designed 
to  compare  the  GA  performance  with  a  random  search.  Case  18  is  an  unconstrained 
case  as  Case  1  (i.e. ,  it  has  the  same  search  space  with  Case  1).  Case  18  was  run  by 
evaluating  (400x600=240000)  random  solutions.  On  the  other  hand,  Case  1  was  run 
with  a  parent  population  number  of  400  and  children  population  number  of  400  for 
600  generations  using  the  GA.  Therefore,  a  maximum  number  of  (400x600=240000) 
solutions  are  evaluated  in  Case  1,  which  is  equal  to  the  number  of  solutions  evaluated 
in  Case  18. 

The  performance  of  Case  1  and  Case  18  are  compared  in  Table  4.8.  Case  18 
has  generated  poor  solutions  at  various  cost  levels.  On  the  other  hand,  the  GA  run 
for  Case  1  has  generated  very  good  and  acceptable  solutions.  In  all  cases,  the  GA 
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resulted  in  more  than  an  order  of  magnitude  in  cost  over  the  purely  random  search. 
Besides,  the  costs  of  the  GA  generated  solutions  are  very  close  to  each  other.  This 
shows  that  the  GA  designed  for  this  study  is  mush  better  than  a  purely  random  search 
tool,  and  it  generates  consistent  solutions.  This  demonstrates  the  value  of  the  GA  for 
efficiently  searching  for  a  solution.  The  generated  solutions  of  Case  1  (GA)  and  Case 
18  (Random)  are  depicted  in  Figure  4.4. 


Table  4.8:  Comparison  of  Case  1  and  Case  18  Runs 


Case  1  (using  GA) 

Case  18  (best  randomly  generated  result) 

Seed 

Av.  GDOP 

Availability(%) 

Cost 

Av.  GDOP 

Availability(%) 

Cost 

1 

4.41 

100 

223.44 

4.81 

61.00 

6141.29 

2 

3.61 

100 

183.79 

7.39 

75.08 

4122.96 

3 

4.45 

99.41 

275.10 

6.75 

73.96 

6279.04 

4 

3.30 

100 

166.71 

8.62 

72.30 

4729.84 

5 

4.05 

100 

206.56 

6.10 

65.36 

5549.88 

4-5  Obtained  Data 

Case  1  is  an  unconstrained  case  and  has  the  largest  search  space.  Having  the 
largest  search  space  among  the  other  cases  ensures  that  the  search  space  of  Case  1 
includes  the  best  possible  solution.  However,  a  larger  search  space  makes  it  harder  to 
find  good  solutions  for  the  GA.  Therefore,  it  shouldn’t  be  expected  that  the  solutions 
generated  for  Case  1  should  be  better  than  the  solutions  generated  for  the  other  cases. 
In  order  to  assist  the  GA  to  find  better  solutions,  the  search  space  was  constrained 
in  various  ways  for  subsequent  cases. 
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Figure  4.4:  Ground  tracks  of  the  generated  solutions  for  Case  1  (1st  column)  and 

Case  18  (2nd  column)  (Random  Seeds  1,2,  and  3  respectively).  No  Constraints. 


First,  a  detailed  analysis  of  the  results  obtained  for  Case  1  is  given.  Next,  the 
obtained  data  for  the  four-satellite  case  is  analyzed.  Finally,  the  rest  of  the  cases  are 
analyzed  in  groups,  which  have  similar  constraints. 

4-5.1  Case  1  (Unconstrained  case).  This  case  was  designed  to  find  a  good 
solution  in  the  whole  search  space  of  this  thesis.  Although  Case  1  has  the  largest 
search  space  in  this  thesis,  the  GA  has  generated  very  good  solutions.  The  GA 
evolution  for  Case  1  was  able  to  converge  in  600  generations.  The  progress  of  the  GA 
evolution  for  Case  1  Seed  4  is  shown  in  Figure  4.6.  The  red  bars  are  the  average  cost 
values  of  the  populations  and  the  black  bars  are  the  cost  value  of  the  best  individual 
in  each  population.  The  progress  of  the  GA  evolution  demonstrates  that  despite  the 
already  accepted  level  of  elitism,  caused  by  the  use  of  the  enlarged  sampling  space,  it 
doesn’t  converge  prematurely. 

The  ground  track  of  the  generated  constellation  for  Case  1  Seed  4,  one  of  the 
best,  is  shown  in  Figure  4.7.  The  average  GDOP  and  the  availability  values  on  top  of 
the  constellation  plots  are  not  weighted  for  the  receivers.  Therefore,  they  shouldn’t  be 
used  for  comparison.  They  are  provided  to  give  an  overall  idea  of  the  performance  of 
the  constellation.  This  condition  applies  to  all  average  GDOP  and  availability  values 
given  in  the  plots  and  tables  except  the  values  for  specific  receiver  sites.  However,  the 
cost  value  can  be  used  for  comparison,  which  is  already  done  to  determine  the  good 
and  bad  constellations  during  the  GA  process. 
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The  geographical  distribution  of  the  average  GDOP  and  availability  are  given 
in  Figure  4.5  for  Case  1  Seed  4  as  an  example.  The  contour  plots  demonstrate  that 
the  performance  of  the  constellation  does  not  have  fluctuations  geographically  in  the 
central  area  of  the  specified  region  (Figure  1.1).  Therefore,  it  can  be  stated  that  the 
constellation  has  a  geographically  smooth  performance. 

The  generated  constellations  for  Case  1  tend  to  exploit  orbits  that  are  combina¬ 
tions  of  low  inclination  orbits,  which  provide  continuous  visibility,  and  high  inclination 
orbits,  which  provide  better  satellite  geometry  (leads  to  lower  GDOP  values).  The 
constellations  generated  for  this  case  provide  100%  availability  and  an  average  GDOP 
value  of  4. 

The  visibility  of  the  satellites  and  the  GDOP  values  over  the  whole  simulation 
time  is  shown  in  Figure  4.8  for  Receiver  Site  1. 

4-5.2  Obtained  Data  from  the  Four- Satellite  Case  (Case  2).  In  theory, 
four  satellites  are  enough  for  creating  a  navigation  satellite  constellation,  as  discussed 
in  Section  2.4.  That  is  still  true,  because  navigation  signals  transmitted  by  four 
satellites  are  generally  sufficient  for  calculating  a  position  estimate.  However,  the 
desired  accuracy  level  may  not  be  achieved  by  a  four-satellite  constellation  if  the 
satellite  geometry  is  poor.  In  order  to  analyze  the  performance  of  an  optimal  four- 
satellite  constellation,  Case  2  was  run  without  any  constraints.  The  GA  runs  for  Case 
2  were  able  to  converge  in  900  generations.  The  evolutionary  progress  of  the  GA  for 
Case  2  Seed  1  is  shown  in  Figure  4.9. 
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Figure  4.5:  a)  Geographical  Distribution  of  Average  GDOP  for  Case  1  Seed  4 
b)  Geographical  Distribution  of  Availability  for  Case  1  Seed  4 
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Figure  4.6:  GA  Evolution  for  Case  1  Seed  4 


Average  GDOP  :  3.30  Availability  :  100.00%  Cost :  166.71 


Figure  4.7:  Ground  Track  of  the  Generated  Constellation  for  Case  1  Seed  4. 
(Average  GDOP  is  not  weighted  for  the  receiver  sites) 
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Satellite  Visibility  at  Site  1 


Mean  GDOP:  3.20  Availability:  100.00  % 


Figure  4.8:  Satellite  Visibility  and  GDOP  at  Receiver  Site  1 
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Figure  4.9:  GA  evolution  of  four-satellite  constellation  (Case  2) 
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The  four-satellite  case  was  run  to  analyze  the  performance  of  a  four-satellite 
navigation  constellation  for  a  mid-latitude  region.  The  obtained  data  yields  that  the 
algorithm  tends  to  generate  constellations  with  a  combination  of  satellite  orbits  with 
very  low  inclination  and  high  inclination  values  in  order  to  have  a  well  distributed 
satellite  geometry.  However,  the  performance  of  the  four-satellite  constellations  are 
very  poor,  which  are  shown  in  Table  4.9. 


Table  4.9:  Summary  of  the  Obtained  Data  for  Case  2 


Case 

Seed 

Cost 

Av.  GDOP 

Availability  (%) 

1 

1814.07 

5.12 

89.88 

2 

2824.25 

6.66 

83.33 

2 

3 

3518.99 

7.99 

79.06 

4 

4067.38 

6.84 

75.51 

5 

4345.97 

5.50 

73.34 

In  the  generated  constellations  for  Case  2,  there  is  continuous  unavailability  for 
around  three  hours  during  the  day.  The  GA-generated  constellations  for  the  four- 
satellite  case  show  that  it  is  hard  to  have  a  good  navigation  performance  with  a 
four-satellite  constellation  over  a  mid-latitude  region.  The  constellations  with  highly 
inclined  satellite  orbits  have  long  periods  of  complete  system  unavailability  when  the 
fourth  satellite  is  not  visible.  The  constellations  with  low  inclined  satellite  orbits  don’t 
have  complete  system  unavailability,  because  the  satellites  are  around  the  GEO  belt 
and  visible.  However,  their  performance  values  are  poor.  The  ground  tracks  of  the 
generated  constellations  for  Case  2  are  shown  in  Figure  4.10. 
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Average  GDOP  :  5.12  Availability  :  89.88%  Cost :  1814.07 


Figure  4.10:  Ground  Tracks  of  the  Generated  Solutions  for  the  Four  Satellite  Con¬ 
stellation  (Case  2) 
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4-5.3  Group  1  ( Case  3,4,5,  and  6).  This  group  of  cases  comprise  of  three 
unconstrained  GSO  satellites  and  two  GEO  satellites  with  different  constraints.  In 
Case  3,  two  GEO  satellites  are  fixed  at  0°  and  70°  degrees  on  the  equator.  In  Case 
4,  two  GEO  satellites  have  twenty  degrees  of  freedom  in  the  range  of  (350:10)  and 
(60:80).  In  Case  5,  the  GEO  satellites  are  free  to  be  placed  anywhere  on  the  equator. 
In  Case  6,  two  GEO  satellites  have  a  freedom  of  120  degrees  on  the  equator  between 
330  and  110  degrees.  The  summary  of  the  obtained  performance  values  for  this  group 
is  given  in  Table  4.10. 


Table  4.10:  Summary  of  the  Obtained  Data  for  Case  3,  4,  5,  and  6 


Case 

Seed 

Cost 

Av.  GDOP 

Availability  (%) 

1 

281.44 

5.55 

100 

2 

264.19 

5.18 

100 

3 

3 

260.73 

5.13 

100 

4 

278.75 

5.48 

100 

5 

275.36 

5.40 

100 

1 

176.24 

3.49 

100 

2 

194.48 

3.81 

100 

4 

3 

178.38 

3.50 

100 

4 

178.09 

3.49 

100 

5 

184.02 

3.61 

100 

1 

196.66 

3.86 

100 

2 

181.35 

3.56 

100 

5 

3 

176.48 

3.46 

100 

4 

161.53 

3.17 

100 

5 

188.25 

3.70 

100 

1 

192.58 

3.78 

100 

2 

193.54 

3.80 

100 

6 

3 

197.17 

3.87 

100 

4 

171.75 

3.37 

100 

5 

231.03 

4.54 

100 
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The  generated  constellation  for  Case  4  Seed  2  looks  like  a  manmade  constellation 
design.  One  of  the  GSO  satellites  was  put  in  a  nearly  circular  and  very  low  inclined 
orbit,  which  is  very  close  to  a  GEO  orbit.  Moreover  this  satellite  was  placed  almost  in 
the  middle  of  the  other  GEO  satellites  on  the  equator.  The  other  GSO  satellites  were 
placed  almost  on  the  same  orbit  with  around  140°  phasing.  The  orbital  parameters 
of  this  particular  constellation  is  given  in  Table  4.11.  The  ground  track  is  depicted  in 
Figure  4.11. 

Table  4.11:  Orbital  Parameters  of  the  Generated  Constellation  for 

Case  4  Seed  2. 


Satellite 

1 

2 

3 

4 

5 

Eccentricity 

0 

0 

0.01 

0.85 

0.85 

Inclination 

0 

0 

0.44 

70.21 

70.12 

LAN 

-10 

80 

33.67 

44.92 

44.95 

Argument  of  Perigee 

0 

0 

90.18 

270.53 

270.53 

True  Anomaly 

0 

0 

360 

319.14 

180.35 

Average  GDOP  :  3.81  Availability  :  100.00%  Cost :  194.48 


Figure  4.11:  Sample  Constellation  of  Case  4  (Seed  2) 
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The  plots  of  the  satellite  visibility  and  GDOP  versus  time  for  this  constellation 
(Case  4  Seed  2)  are  given  in  Figure  4.12  and  4.13.  In  these  plots,  there  are  two 
GDOP  values  3  and  5.  While  there  are  five  visible  satellites,  the  GDOP  value  is 
almost  constant  at  3.  While  there  are  four  visible  satellites,  the  GDOP  value  is 
almost  constant  at  5.  Because  the  three  GEO  satellites  are  stationary  on  the  sky 
and  the  GDOP  value  is  dependent  on  the  position  of  the  GSO  satellites  (cyan(x)  and 
magenta(+)  in  Figure  4.14.d).  This  kind  of  a  constellation  might  be  optimal  for  most 
applications,  because  the  provided  GDOP  value  (which  is  the  only  parameter  that 
determines  the  accuracy  in  this  thesis)  is  quite  predictable.  It  is  approximately  3  or 
5. 

The  common  feature  of  this  group  of  cases  is  two  GEO  satellites  and  three 
unconstrained  GSO  satellites.  In  most  cases,  the  GA  generated  constellations  with 
a  low  inclination  GSO  satellite  and  two  high  inclination  GSO  satellites.  The  low 
inclination  satellite  sometimes  appeared  as  almost  a  GEO  satellite.  This  result  implies 
that  a  constellation  of  three  GEO  satellites  and  two  highly  inclined  GSO  satellites 
might  be  well  suited  for  a  mid-latitude  regional  navigation  satellite  system.  In  Figure 
4.14,  the  ground  tracks  of  the  generated  constellations  with  this  feature  are  given. 

The  constellations  generated  for  this  group  of  cases  usually  exploited  3  GEO 
satellites,  which  provides  continuous  visibility  and  fixed  LOS  vectors,  and  2  GSO 
satellites  with  high  inclination  and  high  eccentricity,  which  provides  good  satellite 
geometry  for  long  periods.  Therefore,  the  acquired  GDOP  values  don’t  have  much 
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Figure  4.12:  Satellite  Visibility  and  GDOP  Plots  for  Case  4  Seed  2  at  Receiver 

Sites  1  through  6  (a  through  f,  respectively). 
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Satellite  Visibility  at  Site  7 


Satellite  Visibility  at  Site  8 
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Figure  4.13:  Satellite  Visibility  and  GDOP 
Sites  7  through  12  (a  through  f,  respectively). 
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Average  GDOP  :  5.55  Availability  :  100.00%  Cost :  281.44 


(a) 


Average  GDOP  :  5.13  Availability  :  100.00%  Cost :  260.73 


(c) 


Average  GDOP  :  3.61  Availability  :  100.00%  Cost:  184.02 
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Average  GDOP  :  3.87  Availability  :  100.00%  Cost :  197.17 


-90 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 

-180  -150  -120  -90  -60  -30  0  30  60  90  120  150  180 

Longitude 


(g) 


Average  GDOP  :  5.19  Availability  :  100.00%  Cost :  264.19 


Longitude 

(b) 


Average  GDOP  :  3.81  Availability  :  100.00%  Cost :  194.48 


Longitude 

(d) 


Average  GDOP  :  3.86  Availability  :  100.00%  Cost :  196.66 


Longitude 

(f) 


Average  GDOP  :  4.54  Availability  :  100.00%  Cost :  231.03 


-60 


-90 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 

-180  -150  -120  -90  -60  -30  0  30  60  90  120  150  180 

Longitude 

(h) 


Figure  4.14:  Ground  Tracks  of  the  Sample  Constellations  with  2  GEO  Satellites 
and  3  Unconstrained  GSO  Satellites. 

a)  Case  3  Seed  1,  b)  Case  3  Seed  2,  c)  Case  3  Seed  3,  d)  Case  4  Seed  2,  e)  Case  4 
Seed  5,  f)  Case  5  Seed  1,  g)  Case  6  Seed  3,  h)  Case  6  Seed  5 
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fluctuation  except  the  time  when  a  satellite  goes  out  sight  and  when  satellite  comes 


into  view. 

4-5.4  Group  2  (Case  7,8,9,  and  10).  In  groups  2,  3,  and  4,  the  inclination 
of  the  GSO  satellites  are  fixed  to  critical  inclination  (63.4°),  which  provides  better 
resistance  against  J2  perturbation  effect  on  the  parameters  of  the  satellite  orbit. 

The  Group  2  cases  are  comprised  of  three  GSO  satellites  with  different  con¬ 
straints  and  two  GEO  satellites  to  be  placed  between  330°  and  110°  on  the  equator. 
In  Case  7,  the  inclination  values  are  fixed  to  63.4°  for  the  GSO  satellites  .  In  Case  8, 
the  inclination  values  and  the  argument  of  perigee  values  are  fixed  to  63.4°  and  270° 
respectively  for  the  GSO  satellites.  In  Case  8,  in  addition  to  the  constraints  of  Case 
7,  the  eccentricity  is  limited  to  be  between  0.3  and  0.7.  In  Case  10,  the  argument  of 
perigee  constraint  is  removed  from  those  of  Case  9.  The  summary  of  the  obtained 
performance  values  for  this  group  is  given  in  Table  4.12. 

In  Case  7,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  85°  separation  on  the  equator.  The  eccentricity  of  the  GSO  satellites  were 
given  values  between  0.6  and  0.8.  The  generated  values  for  the  argument  of  perigee 
of  the  GSO  satellites  were  between  230°  and  315°.  This  range  of  argument  of  perigee 
caused  the  dwell  point  of  the  satellites  to  be  close  to  each  other.  The  ground  track  of 
a  sample  constellation  generated  for  this  case  is  shown  in  Figure  4.15. 

In  Case  8,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  70°  separation  on  the  equator.  The  eccentricity  of  the  GSO  satellites  were 
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Table  4.12:  Summary  of  the  Obtained  Data  for  Case  7,  8,  9,  and  10 


Case 

Seed 

Cost 

Av.  GDOP 

Availability  (%) 

1 

178.34 

3.50 

100 

2 

166.33 

3.26 

100 

7 

3 

168.50 

3.31 

100 

4 

222.65 

4.37 

100 

5 

159.14 

3.13 

100 

1 

247.24 

4.77 

99.97 

2 

193.08 

3.78 

100 

8 

3 

219.52 

4.30 

100 

4 

158.29 

3.62 

100 

5 

190.56 

3.73 

100 

1 

188.01 

3.68 

100 

2 

223.96 

4.39 

100 

9 

3 

181.11 

3.55 

100 

4 

192.45 

3.77 

100 

5 

209.90 

4.11 

100 

1 

158.80 

3.12 

100 

2 

155.22 

3.05 

100 

10 

3 

161.64 

3.17 

100 

4 

161.62 

3.17 

100 

5 

158.92 

3.12 

100 

-180  -150  -120  -90  -60  -30  0  30  60  90  120  150  180 

Longitude 

Figure  4.15:  Ground  Track  of  a  Sample  Constellation  of  Case  7  (Seed  5). 


Average  GDOP  :  3.13  Availability  :  100.00%  Cost :  159.14 
90  i  i  i  i  i  i  i  i  i  i 
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given  values  between  0.55  and  0.85.  The  longitude  of  ascending  nodes  of  the  GSO 
satellites  were  placed  in  an  average  of  20  degree  range.  The  ground  track  of  a  sample 
constellation  generated  for  this  case  is  shown  in  Figure  4.16. 


Average  GDOP  :  3.62  Availability  :  100.00%  Cost :  185.29 


Figure  4.16:  Ground  Track  of  a  Sample  Constellation  of  Case  8  (Seed  4). 


In  Case  9,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  80°  separation  on  the  equator.  The  eccentricity  of  the  GSO  satellites  were 
given  values  between  0.6  and  0.7.  The  longitude  of  ascending  nodes  of  the  GSO 
satellites  were  placed  in  an  average  of  15  degree  range.  The  very  close  longitude  of 
ascending  node  values  for  the  GSO  satellites  almost  put  them  on  the  same  orbit  with 
different  true  anomaly  values.  The  ground  track  of  a  sample  constellation  generated 
for  this  case  is  shown  in  Figure  4.17. 

In  Case  10,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  95°  separation  on  the  equator.  The  eccentricity  of  the  GSO  satellites  were 
given  values  between  0.4  and  0.7.  The  longitude  of  ascending  nodes  of  the  GSO 
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Average  GDOP  :  3.55  Availability  :  100.00%  Cost :  181.11 


Figure  4.17:  Ground  Track  of  a  Sample  Constellation  of  Case  9  (Seed  3). 

satellites  were  placed  in  an  average  of  20  degree  range.  The  argument  of  perigee 
values  of  the  GSO  satellites  were  between  220°  and  270°.  The  ground  track  of  a 
sample  constellation  generated  for  this  case  is  shown  in  Figure  4.18. 


Average  GDOP  :  3.05  Availability  :  100.00%  Cost :  155.22 


Figure  4.18:  Ground  Track  of  a  Sample  Constellation  of  Case  10  (Seed  2). 


In  Group  2,  the  GEO  satellites  are  generated  with  an  average  of  85°  separation 
between  them  on  the  equatorial  plane.  The  GSO  satellites  tended  to  have  similar 
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LAN  and  argument  of  perigee  values.  Moreover,  in  some  constellations,  the  GSO 
satellites  were  placed  in  almost  identical  orbits.  Because  of  the  fixed  inclination 
values  of  the  GSO  satellites,  the  satellite  geometry  became  more  unstable  than  Group 
1  constellations.  However,  the  geometric  diversity  of  three  GSO  satellite  improved 
the  performance. 

4-5.5  Group  3  (Case  11,12,  and  13).  This  group  of  cases  are  comprised 
of  three  GSO  satellites  with  different  constraints  than  Group  2  cases  and  two  GEO 
satellites  to  be  placed  between  330°  and  110°  on  the  equator.  In  Case  11,  the  inclina¬ 
tion  values  and  the  argument  of  perigee  values  are  fixed  to  63.4°  and  270°  respectively 
and  the  longitude  of  ascending  node  values  are  constrained  to  be  between  0°  and  180° 
on  the  equator  for  the  GSO  satellites.  In  Case  12,  unlike  Case  11  the  argument  of 
perigee  values  of  the  GSO  satellites  are  not  fixed  and  variable  between  180°  and  360°. 
In  Case  13,  in  addition  to  the  constraints  of  Case  12,  the  eccentricity  values  of  the 
GSO  satellites  are  limited  to  the  range  of  0.3  to  0.7.  The  summary  of  the  obtained 
performance  values  for  this  group  is  given  in  Table  4.13. 

In  Case  11,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  70°  separation  on  the  equator.  The  eccentricity  of  the  GSO  satellites  were 
given  values  between  0.55  and  0.80.  The  longitude  of  ascending  nodes  of  the  GSO 
satellites  were  placed  in  an  average  of  25  degree  range.  The  very  close  longitude  of 
ascending  node  values  for  the  GSO  satellites  almost  put  them  on  the  same  orbit  with 
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Table  4.13:  Summary  of  the  Obtained  Data  for  Case  11,  12,  and  13 


Case 

Seed 

Cost 

Av.  GDOP 

Availability  (%) 

1 

261.80 

5.13 

100 

2 

272.05 

5.25 

99.97 

11 

3 

227.23 

4.45 

100 

4 

240.16 

4.62 

99.97 

5 

216.03 

4.15 

99.97 

1 

155.82 

3.06 

100 

2 

715.93 

5.32 

97.09 

12 

3 

1276.26 

6.99 

93.83 

4 

288.48 

4.99 

99.76 

5 

1091.96 

6.74 

95.11 

1 

214.77 

3.86 

99.84 

2 

204.92 

4.02 

100 

13 

3 

207.19 

4.06 

100 

4 

239.18 

4.69 

100 

5 

301.47 

5.34 

99.81 

different  true  anomaly  values.  The  ground  track  of  a  sample  constellation  generated 
for  this  case  is  shown  in  Figure  4.19. 


In  Case  12,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  65°  separation  on  the  equator.  The  eccentricity  of  the  GSO  satellites  were 
given  values  between  0.50  and  0.80.  The  argument  of  perigee  values  of  the  GSO 
satellites  were  between  210°  and  300°.  A  sample  ground  track  of  the  constellation 
generated  for  this  case  is  shown  in  Figure  4.20. 


In  Case  13,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  60°  separation  on  the  equator.  The  eccentricity  of  the  GSO  satellites  were 
given  values  between  0.40  and  0.65.  The  argument  of  perigee  values  of  the  GSO 
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Latitude  a  Latitude 


Average  GDOP  :  5.13  Availability  :  100.00%  Cost :  261.80 


Longitude 


4.19:  Ground  Track  of  a  Sample  Constellation  of  Case  11  (Seed  1) 


Average  GDOP  :  3.06  Availability  :  100.00%  Cost :  155.82 


Figure  4.20:  Ground  Track  of  a  Sample  Constellation  of  Case  12  (Seed  1) 


satellites  were  between  230°  and  290°.  A  sample  ground  track  of  the  constellation 


generated  for  this  case  is  shown  in  Figure  4.21. 


Average  GDOP  :  4.02  Availability  :  100.00%  Cost :  204.92 


Figure  4.21:  Ground  Track  of  a  Sample  Constellation  of  Case  13  (Seed  2). 

The  generated  constellations  for  Group  2  performed  well.  The  longitude  of 
ascending  node  values  of  the  generated  GSO  satellites  were  in  the  range  of  0  to  180. 
Therefore,  constraining  the  longitude  of  ascending  node  values  to  this  range  would 
assist  the  GA  to  find  even  better  solutions.  The  cases  of  Group  3  were  designed  with 
this  idea.  However,  the  generated  constellations  for  Group  3  cases  performed  at  the 
same  level  of  Group  2  constellations.  The  best  cost  values  in  the  constellations  of 
Group  2  and  3  are  155.22  (Case  10  Seed  2)  and  155.82  (Case  12  Seed  1).  Those 
constellations  are  also  the  best  of  all  generated  constellations.  Therefore,  the  optimal 
constellation  of  this  problem  might  be  very  similar  those  solutions.  Although  those 
constellations  are  given  in  the  previous  figures,  they  are  given  together  in  Figure  4.22 
with  their  visibility  and  GDOP  versus  time  plots.  For  both  constellations,  the  system 
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is  always  available  ( GDOP  <  10),  and  the  fluctuation  of  the  GDOP  value  within  a 


range  of  2.4  to  3.7.  For  most  applications,  that  much  fluctuation  might  be  within  the 
acceptable  limits. 


Satellite  Visibility  at  Site  1 


Satellite  Visibility  at  Site  1 


Time  (hrs) 

Mean  GDOP:  3.00  Availability:  100.00  % 


Time  (hrs) 

(c) 


Time  (hrs) 

Mean  GDOP:  3.01  Availability:  100.00  % 


0  10  20  30  40  50 
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70  80 


Figure  4.22:  The  Best  Generated  Constellations  in  Group  2  and  3. 

a)  Case  10  Seed  2 

b)  Case  12  Seed  1 


4-5.6  Group  4  (Case  14  and  15).  This  group  of  cases  comprise  of  three 
GSO  satellites  with  some  constraints  and  two  GEO  satellites  to  be  placed  anywhere 
on  the  equator.  In  Case  14,  the  inclination  values  of  the  GSO  satellites  are  fixed  to 
63.4°  and  the  longitude  of  ascending  node  values  are  constrained  to  be  between  0° 
and  180°  on  the  equator  for  the  GSO  satellites.  Moreover,  the  eccentricity  values  and 


the  argument  of  perigee  values  are  constrained  to  be  in  the  range  of  0.3  to  0.7  and 


180°  to  360°  respectively.  In  Case  15,  unlike  Case  14,  the  constraints  for  the  longitude 


of  ascending  node  and  the  argument  of  perigee  values  are  removed.  The  summary  of 


the  obtained  performance  values  for  this  group  is  given  in  Table  4.14. 


Table  4.14:  Summary  of  the  Obtained  Data  for  Case  14  and  15 


Case 

Seed 

Cost 

Av.  GDOP 

Availability  (%) 

1 

259.46 

4.52 

99.79 

2 

167.732 

4.14 

100 

14 

3 

185.17 

3.55 

99.97 

4 

175.15 

3.44 

100 

5 

270.54 

4.82 

99.84 

1 

187.89 

3.68 

100 

2 

186.10 

3.65 

100 

15 

3 

206.61 

4.05 

100 

4 

259.66 

4.77 

99.89 

5 

178.45 

3.50 

100 

In  Case  14,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  80°  separation  on  the  equator.  The  GSO  satellites  are  generated  with 
different  characteristics.  The  ground  track  of  a  sample  constellation  generated  for 
this  case  is  shown  in  Figure  4.23. 

In  Case  15,  the  GEO  satellites  of  the  generated  solutions  were  placed  with  an 
average  of  70°  separation  on  the  equator.  The  GSO  satellites  were  generated  with 
different  characteristics.  The  ground  track  of  a  sample  constellation  generated  for 
this  case  is  shown  in  Figure  4.24. 

This  group  of  cases  were  designed  to  see  if  removing  the  constraint  on  the  longi¬ 
tude  of  ascending  node  values  of  GEO  improves  the  performance  of  the  constellation. 


Average  GDOP  :  3.29  Availability  :  100.00%  Cost :  167.73 


Longitude 

Figure  4.23:  Ground  Track  of  a  Sample  Constellation  of  Case  14  (Seed  2). 


Average  GDOP  :  3.50  Availability  :  100.00%  Cost :  178.45 


Figure  4.24:  Ground  Track  of  a  Sample  Constellation  of  Case  15  (Seed  5). 
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The  separation  between  the  GEO  satellites  of  the  constellations  were  90°  on  average, 
which  is  not  very  different  from  the  values  obtained  in  the  previous  cases.  Besides, 
there  were  no  improvement  in  the  cost  value.  Therefore,  the  optimal  separation  be¬ 
tween  the  GEO  satellites  for  this  problem  might  be  90°. 


4.. 5.1  Group  5  (Case  16  and  17).  This  group  of  cases  comprise  of  five 
GSO  satellites  different  constraints.  I11  Case  16,  only  the  true  anomaly  values  are 
constrained  to  the  range  of  180°  to  360°.  In  Case  17,  in  addition  to  the  constraints  of 
Case  16,  the  inclination  values  and  the  longitude  of  ascending  node  values  are  limited 
to  the  range  of  0°  to  75°  and  0°  to  180°  respectively.  The  summary  of  the  obtained 
performance  values  for  this  group  is  given  in  Table  4.15. 

Table  4.15:  Summary  of  the  Obtained  Data  for  Case  16  and  17 


Case 

Seed 

Cost 

Av.  GDOP 

Availability  (%) 

1 

1652.18 

6.99 

91.56 

2 

2158.40 

11.70 

88.65 

16 

3 

1882.42 

6.20 

89.69 

4 

2222.04 

10.74 

87.77 

5 

2533.44 

5.70 

84.40 

1 

2234.58 

7.60 

86.78 

2 

2487.85 

7.85 

85.63 

17 

3 

685.40 

4.83 

97.14 

4 

2084.21 

6.15 

87.13 

5 

912.81 

7.48 

95.59 

In  Case  16,  various  constellations  were  generated  with  nothing  in  common. 
Their  performance  are  very  poor  than  most  cases.  The  ground  track  of  a  sample 
constellation  for  this  case  is  shown  in  Figure  4.25. 
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Average  GDOP  :  6.99  Availability  :  91 .56%  Cost :  1652.18 


Figure  4.25:  Ground  Track  of  a  Sample  Constellation  of  Case  16  (Seed  1). 

In  Case  17,  the  generated  constellations,  like  Case  16,  had  nothing  in  common 
but  poor  performance.  The  ground  track  of  a  sample  constellation  generated  for  this 
case  is  shown  in  Figure  4.26. 


Average  GDOP  :  4.83  Availability  :  97.14%  Cost :  685.40 


Figure  4.26:  Ground  Track  of  a  Sample  Constellation  of  Case  17  (Seed  3). 

The  cases  of  Group  5  was  designed  to  see  if  the  constellation  is  able  to  demon¬ 
strate  good  performance  without  GEO  satellites.  The  emphasis  given  to  the  com- 
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binations  of  orbits  which  have  the  apogee  point  over  the  northern  hemisphere.  The 
objective  was  to  increase  the  visibility  over  the  receiver  points.  This  constraint  on 
the  argument  of  perigee  shrank  the  possible  volume  of  the  satellite  geometry.  Thus, 
the  generated  constellations  demonstrated  poorer  performance  values  than  the  con¬ 
stellations  of  most  other  cases.  Some  of  the  constellations  of  this  group  performed  at 
the  same  level  of  four-satellite  constellations  (Case  2).  Therefore,  it  is  wise  to  employ 
GEO  satellites  in  a  mid-latitude  regional  navigation  constellation. 

4-5.8  Summary  of  the  Results.  The  summary  of  the  obtained  data  for  all 
cases  are  given  in  tables  4.16,  4.17,  and  4.18.  The  results  show  that  a  four-satellite 
constellation  has  poor  navigation  performance  over  a  mid-latitude  region.  However, 
the  performance  of  the  five-satellite  constellations  are  good.  Therefore,  a  mid-latitude 
regional  navigation  system  should  comprise  of  at  least  five  satellites.  The  cases  with 
GEO  satellites  and  three  unconstrained  GSO  satellites  generated  results  with  three 
GEO  satellites.  Thus,  a  constellation  with  three  GEO  satellites  and  two  GSO  satellites 
appears  to  be  a  reasonable  solution.  Additionally,  the  constellations  which  don’t  have 
GEO  satellites  demonstrated  poor  performance  values.  Therefore,  it  can  be  concluded 
that  the  GEO  satellites  are  indispensable  for  this  problem. 

4- 6  J2  Perturbation  Effect  on  the  Performance  of  the  Constellation 

The  Earth’s  non-spherical  shape  causes  periodic  variations  in  all  of  the  orbital 
elements.  The  dominant  effects  are  in  the  RAAN  and  argument  of  perigee.  The  rate 
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Table  4.16:  The  Summary  of  the  Performance  of  All  Cases  (1). 


Case 

Description 

Cost 

Average 

GDOP 

Availa¬ 

bility 

Bit 

Length 

223.44 

4.41 

100 

183.79 

3.61 

100 

1 

Unconstrained 

275.10 

4.45 

99.41 

235 

166.71 

3.30 

100 

206.56 

4.05 

100 

1814.07 

5.12 

89.88 

Unconstrained 

2824.25 

6.66 

83.33 

2 

3518.99 

7.99 

79.06 

193 

Four-Satellite 

4067.38 

4345.97 

6.84 

5.50 

75.51 

73.34 

281.44 

5.55 

100 

2  GEO  fW(0,70) 

264.19 

5.19 

100 

3 

260.73 

5.13 

100 

151 

3  GSO 

278.75 

5.48 

100 

275.36 

5.40 

100 

176.24 

3.49 

100 

2  GEO  0LAiV(350:10, 60:80) 

194.48 

3.81 

100 

4 

178.38 

3.50 

100 

165 

3  GSO 

178.09 

3.49 

100 

184.02 

3.61 

100 

196.66 

3.86 

100 

2  GEO 

181.35 

3.56 

100 

5 

176.48 

3.46 

100 

173 

3  GSO 

161.53 

3.17 

100 

188.25 

3.70 

100 

192.57 

3.78 

100 

2  GEO  0LAiV(330:110) 

193.54 

3.80 

100 

6 

197.17 

3.87 

100 

171 

3  GSO 

171.75 

3.37 

100 

231.03 

4.54 

100 
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Table  4.17:  The  Summary  of  the  Performance  of  All  Cases  (2). 


Case 

Description 

Cost 

Average 

GDOP 

Availa¬ 

bility 

Bit 

Length 

7 

2  GEO  DLav(330:110) 

3  GSO  7(63.4) 

178.34 

166.33 

168.50 

222.65 

159.14 

3.50 

3.26 

3.31 

4.37 

3.13 

100 

100 

100 

100 

100 

144 

8 

2  GEO  DLAv(330:110) 

3  GSO  f (63.4)  w(270) 

247.24 

193.08 

219.52 

185.29 

190.56 

4.77 

3.78 

4.30 

3.62 

3.73 

99.97 

100 

100 

100 

100 

120 

9 

2  GEO  OLAJV(330:110) 

3  GSO  e(0.3:0.7)  7(63.4)  ca(270) 

188.01 

223.96 

181.11 

192.45 

209.90 

3.68 

4.39 

3.55 

3.77 

4.11 

100 

100 

100 

100 

100 

117 

10 

2  GEO  Dlav(330:110) 

3  GSO  e(0.3:0.7)  7(63.4) 

158.80 

155.22 

161.64 

161.62 

158.92 

3.12 

3.05 

3.17 

3.17 

3.12 

100 

100 

100 

100 

100 

141 

11 

2  GEO  DLAiV(330:110) 

3  GSO  f (63.4)  Dla/v(0:180)  ca(270) 

261.80 

272.05 

227.23 

240.16 

216.03 

5.13 

5.25 

4.45 

4.62 

4.15 

100 

100 

100 

99.97 

99.97 

117 

12 

2  GEO  DLAv(330:110) 

3  GSO  f (63.4)  DiAJV(0:180)  cn(180:360) 

155.82 

715.93 

1276.26 

288.48 

1091.96 

3.06 

5.32 

6.99 

4.99 
6.74 

100 

97.09 

93.83 

99.76 

95.11 

138 

13 

2  GEO  Dlav(330:110) 

3  GSO  e(0.3:0.7)  7(63.4)  DiAJV(0:180) 
c^(180:360) 

214.77 

204.92 

207.19 

239.18 

301.47 

3.86 

4.02 

4.06 

4.69 

5.34 

99.84 

100 

100 

100 

99.81 

135 

95 


Table  4.18:  The  Summary  of  the  Performance  of  All  Cases  (3). 


Case 

Description 

Cost 

Average 

GDOP 

Availa¬ 

bility 

Bit 

Length 

259.46 

4.52 

99.79 

2  GEO 

167.732 

4.14 

100 

14 

185.17 

3.55 

99.97 

137 

3  GSO  e(0.3:0.7)  z(63.4)  DLAJV(0:180) 

175.15 

3.44 

100 

o;(180:360) 

270.54 

4.82 

99.84 

187.89 

3.68 

100 

2  GEO 

186.10 

3.65 

100 

15 

206.61 

4.05 

100 

143 

3  GSO  e(0.3:0.7)  i(63.4) 

259.66 

178.45 

4.77 

3.50 

99.89 

100 

1652.18 

6.99 

91.56 

2158.40 

11.70 

88.65 

16 

5  GSO  cu(180:360) 

1882.42 

6.19 

89.69 

225 

2222.04 

10.74 

87.77 

2533.44 

5.69 

84.40 

2234.58 

7.60 

86.78 

2487.85 

7.85 

85.63 

17 

5  GSO  7(0:75)  DLj4JV(180:360)  w(180:360) 

685.40 

4.83 

97.14 

225 

2084.21 

6.15 

87.13 

912.81 

7.48 

95.59 

6849.91 

4.81 

56.54 

Unconstrained 

4122.94 

7.39 

75.08 

18 

6279.04 

6.75 

60.66 

235 

Random 

4729.84 

8.62 

72.30 

5549.88 

6.10 

65.36 

2  GEO  DLAiV(0,70) 

4881.34 

4773.99 

10.40 

12.97 

71.50 

73.02 

19 

4862.05 

13.53 

72.65 

41 

3  GSO  e(0.5)  z(63.4)  0LA7V(60,80,100) 

4881.34 

10.40 

71.50 

cu(270)  z/(0:360, 0:360, 180) 

4885.03 

10.56 

71.50 

Benchmark 

20 

2  GEO  DLav(0,70) 

3  GSO  e(0.5)  z(63.4)  0LA7V(60,80,100) 

4740.15 

12.09 

73.03 

- 

cu(270)  z/(0:360, 0:360, 180) 
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of  the  changes  are  as  follows  [19]: 


ClJo  =  —  1.5n  J2(;^)2(cosf)(l  —  e2)  2  (4.1) 

a 

ioj2  =  0.75  n  J2(— -)2(4  —  5sin2z)(l  —  e2)-2  (4.2) 

a 

This  perturbation  is  termed  J2  perturbation.  J2  perturbation  disrupts  the  de¬ 
signed  constellation,  if  orbit  maintenance  maneuvers  are  not  made,  which  reduces 
the  usable  lifetime  of  the  satellite.  An  inclination  value  of  90°  removes  the  effect  on 
the  RAAN  (makes  Qj2  zero),  but  this  inclination  is  not  practical  for  geosynchronous 
orbits.  However,  the  effect  on  the  argument  of  perigee  can  be  removed  with  an  incli¬ 
nation  value  which  makes  (4  —  5sin2i)  zero  (also  makes  ujj2  zero).  This  inclination 
value  is  calculated  as  63.4°  or  116.6°  and  termed  as  critical  inclination. 

This  effect  is  not  included  in  the  orbit  propagation  of  the  GA  runs.  However, 
it  would  be  better  to  know  how  robust  the  generated  constellation  against  J2  pertur¬ 
bation  is.  Therefore,  the  analysis  of  the  performance  of  the  generated  constellations 
with  and  without  an  inclination  value  of  63.4°  for  a  30  day  period  was  made.  The 
parameters  of  the  analyzed  constellations  are  given  in  Table  4.19.  The  first  three  and 
last  three  days  of  ground  tracks  of  the  analyzed  constellations  are  depicted  in  Figure 
4.27.  The  satellite  visibility  and  average  GDOP  versus  time  plots  at  Receiver  Site  1 
are  shown  in  Figure  4.28. 

The  analysis  of  the  obtained  data  shows  that  the  constellations  comprising  of 
the  satellites,  which  has  inclination  values  of  63.4°,  are  more  robust  against  J2  per- 


97 


turbation  effect.  The  summary  of  the  J2  effect  on  the  performance  of  the  analyzed 


constellations  are  given  in  Table  4.20. 


Table  4.19:  The  Constellations  for  J2  Analysis 


Constellation  1 

Constellation  2 

Sat 

1 

2 

3 

4 

5 

1 

2 

3 

4 

5 

e 

0 

0 

0.85 

0.05 

0.68 

0 

0 

0.74 

0.63 

0.74 

i 

0 

0 

78.74 

19.62 

60.09 

0 

0 

63.4 

63.4 

63.4 

^LAN 

-10 

80 

58.37 

34.02 

75.25 

-12.49 

74.85 

81.58 

112.44 

150.95 

u) 

0 

0 

247.28 

338.16 

269.82 

0 

0 

315.62 

309.28 

274.76 

V 

0 

0 

7.75 

39.45 

193.748 

0 

0 

357.18 

168.36 

204.31 

Table  4.20:  Summary  of  the  Obtained  Data  for  J2  Effect  Analysis 


Constellation 

Period 

Cost 

Av.  GDOP 

Availability  (%) 

1 

First  3  Days 

177.45 

3.52 

100 

Last  3  Days 

606.67 

4.91 

97.29 

2 

First  3  Days 

166.43 

3.27 

100 

Last  3  Days 

196.91 

3.62 

99.92 

\.l  Summary 

In  this  chapter,  the  obtained  data  and  analysis  of  the  data  are  given.  In  Section 
4.2,  the  settings  and  the  parameters  of  the  GA  were  given.  In  Section  4.4,  the  perfor¬ 
mance  of  the  algorithm  to  find  good  solutions  was  validated.  In  Section  4.3,  a  limited 
version  of  the  search  space  of  the  problem  was  analyzed.  In  Section  4.5,  the  obtained 
data  for  the  groups  of  cases  with  similar  constraints  were  analyzed.  In  Section  4.6, 
the  J2  perturbation  effect  on  the  performance  of  the  generated  constellations  were 


First  3  Day  Average  GDOP  :  3.52  Availability  :  100.00%  Cost :  177.45 
Last  3  Day  Average  GDOP  :  4.91  Availability  :  97.29%  Cost :  606.67 


First  3  Day  Average  GDOP  :  3.27  Availability  :  100.00%  Cost :  166.43 
Last  3  Day  Average  GDOP  :  3.62  Availability  :  99.92%  Cost :  196.91 


(b) 

Figure  4.27:  (a)  Ground  Track  of  First  Constellation. 

(b)  Ground  Track  of  Second  Constellation. 

(Light  colored  ground  tracks  (dots)  show  first  3  days,  dark  colored  ground  tracks 
(stars)  show  last  3  days) 
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(c)  (d) 

Figure  4.28:  (a  and  b)  Satellite  visibility  and  GDOP  plots  of  Constellation  1  for 

the  first  3  and  last  3  days  respectively. 

(c  and  d)  Satellite  visibility  and  GDOP  plots  of  Constellation  2  for  the  first  3  and 
last  3  days  respectively. 
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analyzed.  In  Chapter  V,  the  conclusions  and  the  recommendations  for  future  work  is 
given. 


101 


V.  Conclusions  and  Recommendations 


5. 1  Conclusions 

In  this  thesis  investigation,  the  goal  is  to  design  a  constellation  of  geosyn¬ 
chronous  navigation  satellites  which  maximizes  availability  and  accuracy  over  a  spec¬ 
ified  region  of  the  earth  by  using  the  GA  for  optimization.  In  order  to  accomplish  this 
goal  a  simulation  tool  is  developed  to  find  an  optimum  navigation  satellite  constella¬ 
tion  using  Genetic  Algorithms.  Several  cases  with  different  constraints  are  designed 
and  analyzed  for  navigation  performance  over  a  specified  mid-latitude  region  of  the 
earth.  Most  of  the  generated  results  provide  100%  availability  with  different  average 
GDOP  values.  The  best  generated  constellation  provides  100%  availability  with  an 
average  GDOP  value  of  3.05.  That  performance  level  meets  the  objectives  for  the 
solution  of  this  thesis. 

The  GAs  are  problem  specific  and  need  to  be  designed  according  to  the  prob¬ 
lem.  The  main  objective  of  designing  a  GA  is  to  find  the  global  minimum  of  the 
problem  search  space  while  avoiding  the  local  minimums  in  the  search  space.  That 
can  achieved  by  a  well  designed  GA  which  searches  the  whole  search  space  while  the 
population  evolves.  The  convergence  rate  of  the  GA  evolution  is  a  good  measure  to 
justify  the  quality  of  the  designed  GA.  In  the  executed  runs,  GA  demonstrates  a  good 
convergence  rate. 

The  executed  tests  to  validate  the  use  of  the  GA  to  find  good  solutions  for  this 
problem  shows  that  the  design  and  the  parameters  of  the  GA  are  satisfactory.  The 
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tests  show  that  the  designed  GA  is  far  better  than  a  pure  random  search.  Moreover, 
the  designed  GA  generates  near-optimal  solutions  in  cases  where  the  optimal  solu¬ 
tion  is  known.  Therefore,  the  GA  developed  for  this  study  satisfactorily  meets  the 
objectives. 

The  results  of  the  analysis  of  the  four-satellite  constellation  indicate  that  a 
four-satellite  constellation  is  not  suitable  for  navigation  purposes  over  a  mid-latitude 
region  of  the  earth.  The  generated  solutions  either  had  an  average  of  three-hour 
unavailability  periods  or  very  poor  performance  values,  which  resulted  in  very  high 
cost  values. 

The  cases  with  two  GEO  satellites  and  three  GSO  satellites  without  any  con¬ 
straints  on  the  inclination  values  turned  into  solutions  which  have  three  GEO  satellites 
and  two  highly  inclined  GSO  satellites  after  the  GA  optimization.  The  GEO  satellites 
are  placed  on  the  equator  with  an  average  separation  of  85°,  when  it  is  not  constrained. 
Moreover,  these  kind  of  solutions  performed  very  well  among  the  other  generated  so¬ 
lutions.  Therefore,  it  can  be  stated  that  the  GEO  satellites  are  indispensable  for  this 
problem. 

The  constellations  with  two  GEO  satellites  and  three  GSO  satellites  with  con¬ 
straint  on  the  inclination  value  along  with  the  other  constraints  also  performed  well. 
In  those  cases,  the  GSO  satellites  usually  tended  to  be  phased  on  the  nearly  same 
orbital  plane. 


103 


The  cases  which  are  comprised  of  five  GSO  satellites  are  also  analyzed.  When 
the  apogee  point  of  the  satellite  constrained  to  be  over  the  northern  hemisphere 
(a;  (180  :  360)),  poor  constellations  were  generated.  The  reason  for  the  design  of  cases 
with  this  constraint  was  to  increase  the  visibility  over  the  receiver  points.  However, 
this  constraint  shrank  the  limits  of  possible  satellite  geometry  which  eventually  led 
to  high  GDOP  values. 

The  analysis  of  the  effect  of  J2  perturbation  was  also  made.  The  results  of 
this  analysis  show  that  the  constellations  comprising  of  GEO  satellites  and  critically 
inclined  (i  =  63.4)  GSO  satellites  are  more  robust  against  J2  perturbation.  Besides, 
those  constellations  have  good  navigation  performance  values. 

5. 2  Recommendations 

The  design  of  the  cost  function  is  vital  for  the  performance  of  the  GA.  Because 
the  cost  function  is  used  to  distinguish  the  good  solutions  among  the  other  generated 
solutions,  a  badly  designed  cost  function  may  prevent  the  GA  from  finding  good 
solutions.  Therefore,  the  cost  function  should  state  the  exact  performance  criteria 
of  the  problem.  In  this  research,  the  designed  cost  function  performed  well  in  terms 
finding  good  results.  Because  the  criteria  applied  to  the  cost  function  are  average 
GDOP  value  and  average  availability,  the  algorithm  generated  availability  as  high  as 
possible  and  average  GDOP  as  low  as  possible  values.  However,  for  some  applications, 
the  satellite  navigation  system  should  provide  stable  GDOP  values  along  with  low 
GDOP  average  values.  Therefore,  an  added  criteria  in  the  cost  function  might  be  to 
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minimize  the  magnitude  of  the  peaks  of  GDOP  values.  That  kind  of  cost  function 
may  generate  solution  with  more  stable  performance. 

As  it  was  discussed  in  Section  3.3,  enlarged  sampling  space  is  used  in  the  de¬ 
signed  GA.  Therefore,  a  level  of  elitism  was  already  accepted,  which  might  have  re¬ 
sulted  early  convergence.  In  order  to  mitigate  this  effect  of  enlarged  sampling  space, 
fewer  parent  population  size  and  more  children  population  size  might  be  used  for  this 
problem. 
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Appendix  A.  Ground  Tracks  of  the  Generated  Constellations 
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Average  GDOP  :  4.05  Availability  :  100.00%  Cost :  206.56 
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Figure  A.l:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  1 

(Seed  1  trough  5,  respectively) 
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Average  GDOP  :  5.12  Availability  :  89.88%  Cost :  1814.07 


(a) 


Average  GDOP  :  7.99  Availability  :  79.06%  Cost :  351 8.99 
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Average  GDOP  :  6.66  Availability  :  83.33%  Cost :  2824.25 
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Figure  A. 2:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  2 

(Seed  1  through  5,  respectively) 
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Average  GDOP  :  5.55  Availability  :  100.00%  Cost :  281.44 
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Average  GDOP  :  5.13  Availability  :  100.00%  Cost :  260.73 
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Average  GDOP  :  5.48  Availability  :  1 00.00%  Cost :  278.75 
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Figure  A. 3:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  3 

(Seed  1  through  5,  respectively) 
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Average  GDOP  :  3.49  Availability  :  100.00%  Cost :  176.24 
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Average  GDOP  :  3.50  Availability  :  100.00%  Cost :  178.38 
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Figure  A. 4:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  4 

(Seed  1  through  5,  respectively) 
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Average  GDOP  :  3.86  Availability  :  100.00%  Cost :  196.66 
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Figure  A. 5:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  5 

(Seed  1  through  5,  respectively) 
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Average  GDOP  :  3.78  Availability  :  100.00%  Cost :  192.57  Average  GDOP  :  3.80  Availability  :  100.00%  Cost :  193.54 
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Figure  A. 6:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  6 

(Seed  1  through  5,  respectively) 
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Latitude  Latitude 


Average  GDOP  :  3.50  Availability  :  100.00%  Cost :  178.34 


Average  GDOP  :  3.26  Availability  :  1 00.00%  Cost :  1 66.33 
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Figure  A. 7:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  7 

(Seed  1  through  5,  respectively) 
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Latitude  Latitude 


Average  GDOP  :  4.77  Availability  :  99.97%  Cost :  247.24  Average  GDOP  :  3.78  Availability  :  100.00%  Cost :  193.08 
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Figure  A. 8:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  8 

(Seed  1  through  5,  respectively) 
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Latitude  Latitude 


Average  GDOP  :  3.68  Availability  :  100.00%  Cost :  188.01 
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Average  GDOP  :  3.55  Availability  :  100.00%  Cost :  181.11 
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(c) 


Average  GDOP  :  4.39  Availability  :  1 00.00%  Cost :  223.96 
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Figure  A. 9:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  9 

(Seed  1  through  5,  respectively) 
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Latitude  Latitude 


Average  GDOP  :  3.12  Availability  :  100.00%  Cost :  158.80 
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Figure  A.  10:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  10 
(Seed  1  through  5,  respectively) 
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Average  GDOP  :  4.15  Availability  :  99.97%  Cost :  216.03 
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Figure  A. 11:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  11 
(Seed  1  through  5,  respectively) 
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Figure  A.  12:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  12 
(Seed  1  through  5,  respectively) 
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Average  GDOP  :  3.86  Availability  :  99.84%  Cost:  214.77 
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Figure  A.  13:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  13 
(Seed  1  through  5,  respectively) 
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Average  GDOP  :  4.52  Availability  :  99.79%  Cost :  259.46 
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Figure  A.  14:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  14 
(Seed  1  through  5,  respectively) 


119 


Latitude  Latitude 


Average  GDOP  :  3.68  Availability  :  100.00%  Cost :  187.89  Average  GDOP  :  3.65  Availability  :  100.00%  Cost :  186.10 
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Figure  A.  15:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  15 
(Seed  1  through  5,  respectively) 
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Average  GDOP  :  6.99  Availability  :  91.56%  Cost :  1652.18  Average  GDOP  :  11.70  Availability  :  88.65%  Cost :  2158.40 
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Figure  A.  16:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  16 
(Seed  1  through  5,  respectively) 
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Average  GDOP  :  7.60  Availability  :  86.78%  Cost :  2234.58  Average  GDOP  :  7.85  Availability  :  85.63%  Cost :  2487.85 
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Figure  A.  17:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  17 
(Seed  1  through  5,  respectively) 
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Latitude  Latitude 


Average  GDOP  :  4.81  Availability  :  56.54%  Cost :  6849.91  Average  GDOP  :  7.39  Availability  :  75.08%  Cost:  4122.94 
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Figure  A.  18:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  18 
(Seed  1  through  5,  respectively) 
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Average  GDOP  :  13.53  Availability  :  72.65%  Cost :  4862.05 
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Figure  A.  19:  a  through  e)  Ground  Tracks  of  the  Generated  Solutions  for  Case  19 
(Seed  1  through  5,  respectively) 
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Figure  A. 20: 


Average  GDOP  :  12.09  Availability  :  73.02%  Cost :  4740.15 


Ground  Tracks  of  the  Generated  Solutions  for  Case  20  (Benchmark) 
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